# Comparison of Galerkin, POD and Nonlinear-Normal-Modes Models for Nonlinear Vibrations of Circular Cylindrical Shells

Marco Amabili, Cyril Touzé, Olivier Thomas

► To cite this version:

Marco Amabili, Cyril Touzé, Olivier Thomas. Comparison of Galerkin, POD and Nonlinear-Normal-Modes Models for Nonlinear Vibrations of Circular Cylindrical Shells. ASME IMECE 2006, International Mechanical Engineering Congress and Exposition, 2006, Chicago, United States. pp.373-382, 10.1115/IMECE2006-14602. hal-01154712

HAL Id: hal-01154712

https://ensta-paris.hal.science/hal-01154712

Submitted on 27 May 2015

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
COMPARISON OF GALERKIN, POD AND NONLINEAR-NORMAL-MODES MODELS
FOR NONLINEAR VIBRATIONS OF CIRCULAR CYLINDRICAL SHELLS

M. Amabili

Dipartimento di Ingegneria Industriale,
Università' di Parma, Area delle Scienze 181/a,
Parma, 43100 Italy
marco@me.unipr.it

C. Touzé

ENSTA-UME, Unité de Recherche en Mécanique,
Chemin de la Hunière,
91761 Palaiseau Cedex, France
Cyril.Touze@ensta.fr

O. Thomas

Structural Mechanics and Coupled Systems Laboratory,
CNAM, 2 rue Conté,
75003 Paris, France
olivier.thomas@cnam.fr

ABSTRACT

The aim of the present paper is to compare two different
methods available to reduce the complicated dynamics
exhibited by large amplitude, geometrically nonlinear
vibrations of a thin shell. The two methods are: the proper
orthogonal decomposition (POD) and an asymptotic
approximation of the Nonlinear Normal Modes (NNMs) of the
system. The structure used to perform comparisons is a water-
filled, simply supported circular cylindrical shell subjected to
harmonic excitation in the spectral neighbourhood of the
fundamental natural frequency. A reference solution is obtained
by discretizing the Partial Differential Equations (PDEs) of
motion with a Galerkin expansion containing 16 eigenmodes.
The POD model is built by using responses computed with the
Galerkin model; the NNM model is built by using the
discretized equations of motion obtained with the Galerkin
method, and taking into account also the transformation of
damping terms. Both the POD and NNMs allow to reduce
significantly the dimension of the original Galerkin model. The
computed nonlinear responses are compared in order to verify
the accuracy and the limits of these two methods. For vibration
amplitudes equal to 1.5 times the shell thickness, the two
methods give very close results to the original Galerkin model.
By increasing the excitation and vibration amplitude,
significant differences are observed and discussed.

Keywords: nonlinear vibrations, POD, nonlinear normal modes, shells.

1. INTRODUCTION

Reduced-order models are an attractive research topic in
nonlinear dynamics of fluid and solid systems. From far, the
two most popular methods used to build reduced-order models
(ROMs) are the proper orthogonal decomposition (POD) and
the nonlinear normal modes (NNMs) methods. The first one
(POD, also referred to as the Karhunen-Loève method) uses a
cloud of points in phase space, obtained from simulations or
from experiments, in order to build the reduced subspace that
will contain most information (Zahorian and Rothenberg 1981;
Aubry et al. 1988; Sirovich 1987; Breuer and Sirovich 1991;
Azeez and Vakakis 2001; Sarkar and Païdoussis 2003;
Kerschen, et al. 2003; Amabili et al. 2003; Sarkar and
Païdoussis 2004; Kerschen, et al. 2005; Amabili et al. 2006).
The method is, in essence, linear, as it furnishes the best
orthogonal basis, which decorrelates the signal components and
maximizes variance.

Amabili et al. (2003, 2006) compared Galerkin and POD
models of a water-filled circular cylindrical shell from
moderate to extremely large vibration amplitudes. Accurate
POD models can be built by using only POD modes with
significant energy. In particular, Amabili et al. (2006) found
that more proper orthogonal modes are necessary to reach
energy convergence using time series extracted from more
complex responses (chaotic or quasi-periodic) than from the
periodic ones. Therefore by using complex response it is
possible to build models with larger dimension, suitable to
describe with accuracy large variations of the system parameters.

The second method (NNMs), allows to construct and define the researched subspaces from specific properties of the dynamical systems, by adapting the reduction theorems provided by the mathematics: center manifold theorem (Carr 1981; Guckenheimer and Holmes 1983) and normal form theory (Poincaré 1892; Iooss and Adelmeyer 1998; Elphick et al. 1987). Their application to vibratory systems led to two definitions of Nonlinear Normal Modes, which are equivalent in a conservative framework: either a family of periodic orbits in the vicinity of the equilibrium point (Rosenberg 1966; Mikhlin 1995; Vakakis et al. 1996), or an invariant manifold containing these periodic orbits (Shaw and Pierre, 1991). Numerous asymptotic methods have been proposed for their computation, by application of the center manifold theorem (Shaw and Pierre 1993), the normal form theory (Jézéquel and Lamarque, 1991, Touzé et al. 2004, Touzé and Amabili 2005), the conservation of energy for conservative systems (King and Vakakis 1994), or the method of multiple scales (Lacarbonara et al., 2003). Numerical procedures have also been proposed, recently by Jiang et al (2005a, b), who extended the method proposed by Pescheck et al. (2002) for conservative cases. Bellizzi and Bouc (2005) propose a numerical resolution of an extended KBM (Krylov-Bogoliubov-Mitropolsky) method, while Slater (1996) used continuation techniques to generate the NNM.

Application of the POD and NNMs methods to reduced-order modelling enabled to show that a few degrees of freedom are generally enough to catch the nonlinear behaviour of many structures, versus the several necessary in the corresponding Galerkin models.

The aim of the present study is to provide a full comparison of the results given by the two reduction methods on a realistic example: a water-filled circular cylindrical shell. The reference solution is obtained by the Galerkin method. Its convergence has been carefully checked (Pellicano et al. 2002), and experimental comparisons have been performed (Amabili 2003). The construction of the POD model has been exhaustively explained in Amabili et al. (2003, 2006), whereas the asymptotic NNMs procedure used here is fully explained in Touzé and Amabili (2005). The peculiarity of the NNMs formulation is that damping is taken into account via an improvement of the real normal form calculation presented in Touzé et al. (2004). Comparisons will be drawn on two different cases. First, the ability of the methods to recover frequency-response curves will be investigated, for moderate values of the amplitude of the external force. Then, bifurcation diagrams for varying amplitude of the forcing, leading to complex dynamics, will be discussed.

## 2. BASIC EQUATIONS FOR NONLINEAR VIBRATIONS OF FLUID-FILLED SHELLS

A cylindrical coordinate system (O; x, r, θ) is chosen, with the origin O placed at the centre of one end of the shell. The displacements of the middle surface of the shell are denoted by

u, v and w, in the axial, circumferential and radial directions, respectively; w is taken positive inwards. By using Donnell's nonlinear shallow-shell theory, the equation of motion for finite-amplitude transverse deflection is given by Amabili and Païdoussis (2003)

$$D\nabla^4 w + c h \dot{w} + \rho h \ddot{w} = f - p + \frac{1}{R} \frac{\partial^2 F}{\partial x^2} + \frac{1}{R^2} \left[ \frac{\partial^2 F}{\partial \theta^2} \left( \frac{\partial^2 w}{\partial x^2} \right) - 2 \frac{\partial^2 F}{\partial x \partial \theta} \left( \frac{\partial^2 w}{\partial x \partial \theta} \right) + \frac{\partial^2 F}{\partial x^2} \left( \frac{\partial^2 w}{\partial \theta^2} \right) \right], \quad (1)$$

where $D = E h^3 / [12(1-ν^2)]$ is the flexural rigidity, E Young's modulus, ν the Poisson ratio, h the shell thickness, R the mean shell radius, ρ the mass density of the shell, c the coefficient of viscous damping, p the radial pressure applied to the surface of the shell by the contained fluid, and f is an external local excitation:

$$f = \tilde{f} \delta(R\theta - R\tilde{\theta}) \delta(x - \tilde{x}) \cos(\omega t), \quad (2)$$

where δ is the Dirac delta function, $\tilde{f}$ is the magnitude of the localized (point) force, and $\tilde{\theta}$ and $\tilde{x}$ give the angular and axial coordinates of the point of application of the force, respectively. The viscous damping model introduced in equation (1) is unrealistic; it will later be replaced by modal damping coefficients in the equations of motion.

In equation (1) the overdot denotes a time derivative, and F is the in-plane Airy stress function, which is given by the following compatibility equation (Amabili and Païdoussis 2003):

$$\frac{1}{Eh} \nabla^4 F = -\frac{1}{R} \frac{\partial^2 w}{\partial x^2} + \frac{1}{R^2} \left[ \left( \frac{\partial^2 w}{\partial x \partial \theta} \right)^2 - \frac{\partial^2 w}{\partial x^2} \frac{\partial^2 w}{\partial \theta^2} \right]. \quad (3)$$

In equations (1) and (3), the biharmonic operator is defined as
$$\nabla^4 = \left[ \frac{\partial^2}{\partial x^2} + \frac{\partial^2}{(R^2 \partial \theta^2)} \right]^2 .$$
Donnell's nonlinear shallow-shell equations are accurate only for modes with $n \ge 4$. Attention is focused on simply supported, circumferentially closed circular cylindrical shells of length L. Moreover, u, v and w must be continuous in θ.

Excitations with frequency close to the natural frequency of the lowest modes of the shell are considered; low-frequency modes are associated with predominantly radial motion and are identified by the pair (m, n), where m is the number of axial half-waves and n is the number of circumferential waves.

The contained fluid is assumed to be incompressible, inviscid and irrotational, so that potential theory can be used to describe fluid motion. Liquid-filled shells vibrating in the low-frequency range satisfy the incompressibility hypothesis very well. Nonlinear effects in the dynamic pressure and in the boundary conditions at the fluid-structure interface are neglected. The shell prestress due to the fluid weight is also neglected. The fluid motion is described by the velocity potential $\Phi$, which satisfies the Laplace equation; cavitation is assumed not to occur. Both ends of the fluid volume,
corresponding to the shell edges, are assumed to be open, so that a zero pressure is imposed there; this physically corresponds to a long shell periodically supported (e.g. with ring stiffeners) or it approximates a shell closed by very thin circular plates. The dynamic pressure p exerted by the contained fluid on the shell is given by Amabili (2003)

$$p = \rho_F (\dot{\Phi})_{r=R} = \sum_{m=1}^{M} \sum_{n=1}^{N} \rho_F [\ddot{A}_{m,n}(t)\cos(n\theta) + \ddot{B}_{m,n}(t)\sin(n\theta)] \times \frac{I_n(\lambda_m R)}{\lambda_m I'_n(\lambda_m R)} \sin(\lambda_m x) \quad (4)$$

where $\rho_F$ is the mass density of the internal fluid, $I_n$ is the modified Bessel function of order n, and $I'_n$ its derivative with respect to the argument.

# 3. REFERENCE SOLUTION AND REDUCED-ORDER MODELS

## 3.1 GALERKIN METHOD

The Galerkin method, employing any set of basis functions $\phi_i$, approximates the nonlinear partial differential equation (PDE) by transforming it into a finite set of coupled ordinary differential equations (ODEs), with the solution being expressed as

$$w(\xi, t) = \sum_{i=1}^{K} q_i(t) \phi_i(\xi), \quad (5)$$

where $t$ is time, $\xi$ is the vector of spatial coordinate $(x, \theta)$ describing the shell middle surface $\Omega$, $q_i(t)$ are the generalized coordinates, and $K$ is the number of generalized coordinates (degrees of freedom), i.e. the number of basis functions assumed. The linear modal base is the best choice for discretizing the shell, as these are the eigenfunctions of the linear operator of the PDE. The orthogonality property of the eigenmodes allows decoupling the ODEs at the linear stage. Other sets of basis functions may be used, with the consequence that the ODEs are linearly coupled, and more functions are needed to attain convergence. The key question in the Galerkin method is the convergence of the solution. In order to have a reasonable number of degrees of freedom, it is important to use the most significant modes. In addition to the asymmetric mode directly driven into vibration by the excitation (driven mode), it is necessary to consider (i) the orthogonal mode having the same shape and natural frequency but rotated by $\pi/(2n)$ (companion mode), (ii) additional asymmetric modes, and (iii) axisymmetric modes. In fact, it has clearly been established that, for large-amplitude shell vibrations, the deformation of the shell involves significant axisymmetric oscillations inwards. According to these considerations, the radial displacement $w$ is expanded by using the eigenmodes of the empty shell, which are unchanged for the completely filled shell (Amabili 2003):

$$w(x, \theta, t) = \sum_{m=1}^{3} \sum_{k=1}^{3} \left[ A_{m,kn}(t) \cos(kn\theta) + B_{m,kn}(t) \sin(kn\theta) \right] \times \sin(\lambda_m x) + \sum_{m=1}^{4} A_{(2m-1),0}(t) \sin(\lambda_{(2m-1)} x), \quad (6)$$

where $kn$ is the number of circumferential waves, $m$ is the number of longitudinal half-waves (only odd values are used for symmetry), $\lambda_m = m\pi/L$, and $t$ is the time; $A_{m,n}(t)$, $B_{m,n}(t)$ and $A_{m,0}(t)$ are the generalized coordinates that are unknown functions of $t$; the mode driven in resonance is $(1, n)$, i.e. the mode for $m = k = 1$. The number of degrees of freedom used in the present numerical calculations is 16 (Amabili 2003).

The presence of pairs of modes having the same shape but different angular orientations, the first one described by $\cos(n\theta)$ (driven mode for the excitation given by equation (2)) and the other by $\sin(n\theta)$ (companion mode), in the periodic response of the shell leads to the appearance of a travelling wave around the shell in the $\theta$ direction when both these modes are active and when they have a relative time shift. This phenomenon is due to the axial symmetry of the system.

When the excitation has a frequency close to the resonance of a particular mode, say $(m = 1, n)$, results for relatively low amplitude excitation (case of periodic response) show that (i) the generalized coordinates $A_{1,n}(t)$ and $B_{1,n}(t)$ have the same frequency as the excitation, (ii) the coordinates $A_{1,2n}(t)$, $B_{1,2n}(t)$, $A_{3,2n}(t)$, $B_{3,2n}(t)$ and all the coordinates associated with axisymmetric modes have twice the frequency of the excitation, and (iii) the coordinates $A_{3,n}(t)$, $B_{3,n}(t)$, $A_{1,3n}(t)$, $B_{1,3n}(t)$, $A_{3,3n}(t)$ and $B_{3,3n}(t)$ have three times the frequency of the excitation.

Expansion (6) used for the radial displacement $w$ satisfies identically the out-of-plane boundary conditions and the continuity of the circumferential displacement. The boundary conditions for the in-plane displacements give very complex expressions when transformed into equations involving $w$. Therefore, they are modified into simpler integral expressions that satisfy boundary conditions on the average (Amabili 2003).

When the expansion of $w$, equations (6), is substituted in the right-hand side of equation (3), a partial differential equation for the stress function $F$ is obtained, composed of the homogeneous and the particular solution.

Equations (4) and equation (6) present the same spatial distribution on the shell surface. Therefore, fluid pressure gives only an inertial effect, which is different for each mode of the expansion. Hence, the fluid is expected to change the nonlinear behaviour of the fluid-filled shell, as a consequence of the fundamental interaction among asymmetric and the axisymmetric modes. Usually the inertial effect of the fluid is larger for axisymmetric modes, thus enhancing the nonlinear softening type behaviour of the shell.

By use of the Galerkin method, 16 second-order, ordinary, coupled nonlinear differential equations are obtained for the variables $A_{m,kn}(t)$, $B_{m,kn}(t)$ and $A_{m,0}(t)$, for $m = 1, \dots, M$ and $k = 1, \dots, 3$, by successively weighting the original equation (1) with
the functions that describe the shape of the modes retained in equation (6). These equations have very long expressions containing quadratic and cubic nonlinear terms and have been obtained by using the Mathematica 4 computer software (Wolfram 1999), in order to perform analytical integrals of trigonometric functions. The generic jth Lagrange equation is divided by the modal mass associated with $\ddot{q}_j$, taking the following form

$$ \ddot{q}_j + \omega_j^2 q_j + 2\zeta_j \omega_j \dot{q}_j + \sum_{i=1}^{K} \sum_{p=1}^{K} g_{ip}^{(j)} q_i q_p + \sum_{i=1}^{K} \sum_{p=1}^{K} \sum_{k=1}^{K} h_{ipk}^{(j)} q_i q_p q_k = \hat{f}_j \cos(\omega t), \quad \text{for } j = 1, \dots, K, \quad (7) $$

where $\hat{f}_j$ is projection of the nondimensionalized force, which must be set equal to zero in all the equations where $q_j = B_{m,k,n}$, $g_{ip}^{(j)}$ are coefficients of quadratic terms and $h_{ipk}^{(j)}$ are coefficients of cubic terms; $\zeta_j$ is the modal damping ratio, substituting here the unrealistic viscous damping introduced in equation (1). In equation (7) each generalized coordinate $q_j$ (and therefore the modal damping $\zeta_j$ and circular frequency

$\omega_j$) has to be referred to mode $(m, n)$, i.e. $q_j = A_{m,n}$ or $B_{m,n}$. For computational convenience a non-dimensionalization of variables is also performed: the time is divided by the period of the resonant mode and the vibration amplitudes are divided by the shell thickness $h$. It can be observed that nonlinear terms do not involve time derivative of $q_j$. By introducing a dummy variable, the K second-order equations are transformed into 2×K first-order nonlinear differential equations that are studied by using (i) the software AUTO 97 (Doedel et al. 1998) for continuation and bifurcation analysis of nonlinear ordinary differential equations, and (ii) direct integration of the equations of motion by using the DIVPAG routine of the Fortran library IMSL. Continuation methods allow following the solution path, with the advantage that unstable solutions can also be obtained; these are not ordinarily attainable by using direct numerical integration. The software AUTO 97 is capable of continuation of the solution, bifurcation analysis and branch switching by using pseudo-arclength continuation and collocation methods.

## 3.2 PROPER ORTHOGONAL DECOMPOSITION (POD) METHOD

The POD method optimally extracts the spatial information necessary to characterize the spatio-temporal complexity and inherent dimension of a system, from a set of temporal snapshots of the response, gathered from either numerical simulations or experimental data. In the present context, the temporal responses are obtained via the conventional Galerkin solution. The proper orthogonal modes obtained by the POD method will be used as a basis in conjunction with the Galerkin

approach. The solution can be expressed by using the basis of the proper orthogonal modes $\psi_i(\xi)$,

$$ w(\xi, t) = \sum_{i=1}^{\tilde{K}} a_i(t) \psi_i(\xi), \quad (8) $$

where $a_i$ are the proper orthogonal coordinates and $\tilde{K}$ is the number of proper orthogonal modes (degrees of freedom) used to build the POD model (in general, significantly lower than K in equation (5) necessary for the conventional Galerkin method).

The displacement field $w$ is divided into its time-mean value $\bar{w}(\xi)$ and the zero-mean response $\tilde{w}(\xi, t) = (w(\xi, t) - \bar{w}(\xi))$. In the POD method, the proper orthogonal modes are obtained by minimizing the objective function

$$ \tilde{\lambda} = \left\langle \left( \psi(\xi) - \tilde{w}(\xi, t) \right)^2 \right\rangle \quad \forall \xi \in \Omega, \quad (9) $$

with $\langle \rangle$ denoting the time-averaging operation and $\psi(\xi)$ the generic POD mode. Minimizing of the objective function (9) is obtained, after some mathematics, by solving the following eigenvalue problem:

$$ \int_{\Omega} \langle \tilde{w}(\xi, t) \tilde{w}(\xi', t) \rangle \psi(\xi') d\xi' = \lambda \psi(\xi), \quad (10) $$

where $\langle \tilde{w}(\xi, t) \tilde{w}(\xi', t) \rangle$ is the time-averaged spatial autocorrelation function.

A Galerkin projection scheme for determining proper orthogonal modes (POMs) semi-analytically, and in parallel to approximate the solution of the PDE has been developed in (Sarkar and Païdoussis 2003).

The optimal number of terms $\tilde{K}$ to be retained can be estimated by $\sum_{i=1}^{\tilde{K}} \lambda_i / \sum_{i=1}^{K} \lambda_i \ge 0.99$ in equation (10); however, for each problem this cut-off value can be different. It would be useful to check the convergence of the solution by increasing the value $\tilde{K}$; over a certain value, the results can become less accurate, because the additional terms introduced in the expansion may be highly noise-polluted.

In some applications, it may be better to use time responses obtained for different system parameters in order to produce better proper orthogonal modes.

The POD method gives a set of equations of motion with exactly the same structure of those obtained with the conventional Galerkin method, equation (7), with the difference of a reduced order of degrees of freedom, i.e. $\tilde{K} < K$. Both travelling waves around the shell in opposite directions have been taken into account (one obtained by direct integration and the other by changing the sign to the generalized coordinates associated to $\sin(n\theta)$ terms in equation (6)) to construct the POD reduced-order model; this is fundamental to reproduce the axisymmetry of the shell with no preferential direction.
### 3.3 ASYMPTOTIC NONLINEAR NORMAL MODES (NNMS) METHOD

Nonlinear Normal Modes (NNMs) are here defined as invariant manifold of the state space. They are moreover chosen tangent at the origin, which corresponds to the position of the structure at rest. An asymptotic procedure, based on the normal form theory, is used to compute the NNMs of the system. The method is here briefly recalled, the interested reader is referred to Touzé and Amabili (2005) for a complete description. In particular, it allows to take modal damping into account in the derivations, hence extending previous results obtained for conservative systems (Touzé et al. 2004). A third-order asymptotic development is applied, in order to perform a nonlinear change of coordinates for the system of damped unforced equation of motion, corresponding to equation (7) with right-hand-side equal to zero. A real formulation is used, so that normal forms are expressed with oscillators. The dummy variable $y_j = \dot{q}_j$ for the nondimensional velocity permits to recast the system of equations into first-order. The nonlinear change of coordinates is:

$$ 
\begin{aligned}
q_j = r_j &+ \sum_{i=1}^{K} \sum_{p \ge i} (a_{ip}^{(j)} r_i r_p + b_{ip}^{(j)} s_i s_p) \\
&+ \sum_{i=1}^{K} \sum_{p=1}^{K} c_{ip}^{(j)} r_i s_p + \sum_{i=1}^{K} \sum_{p \ge i} \sum_{k \ge p} (d_{ipk}^{(j)} r_i r_p r_k + e_{ipk}^{(j)} s_i s_p s_k) \\
&+ \sum_{i=1}^{K} \sum_{p=1}^{K} \sum_{k \ge p} (t_{ipk}^{(j)} s_i r_p + u_{ipk}^{(j)} r_i s_p),
\end{aligned}
$$

$$ 
\begin{aligned}
y_j = s_j &+ \sum_{i=1}^{K} \sum_{p \ge i} (\alpha_{ip}^{(j)} r_i r_p + \beta_{ip}^{(j)} s_i s_p) \\
&+ \sum_{i=1}^{K} \sum_{p=1}^{K} \gamma_{ip}^{(j)} r_i s_p + \sum_{i=1}^{K} \sum_{p \ge i} \sum_{k \ge p} (\lambda_{ipk}^{(j)} r_i r_p r_k + \mu_{ipk}^{(j)} s_i s_p s_k) \\
&+ \sum_{i=1}^{K} \sum_{p=1}^{K} \sum_{k \ge p} (v_{ipk}^{(j)} s_i r_p r_k + \xi_{ipk}^{(j)} r_i s_p s_k),
\end{aligned}
\quad (11) $$

where $r_j$ is the transformed nondimensional displacement and $s_j$ is the transformed nondimensional velocity; other symbols are the transformation coefficients, known by analytical expressions. After substitution of (11) into (7), the dynamics, written with the newly introduced variables ($r_j, s_j$), is expressed in an invariant-based span of the state space by cancelling all the invariant-breaking terms (i.e. the non-resonant terms). As a result, proper truncations can now be realized, as all invariant-breaking terms between oscillators in Eqs. (7), have been cancelled. The reduction can now be applied by simply selecting the most important normal coordinates for simulation (master coordinates), and cancelling all the others. In the case considered here, the minimum model must retain the NNMs corresponding to the driven mode ($r_1, s_1$, that are the continuations of $A_{1,5}$ and $\dot{A}_{1,5}$) and the companion mode ($r_2, s_2$,

that are the continuations of $B_{1,5}$ and $\dot{B}_{1,5}$), as these two modes have the same eigenfrequency (1:1 internal resonance). Finally, the reduced-order model (ROM) built by selecting these two pairs of coordinates writes:

$$ 
\begin{aligned}
\ddot{r}_1 + \omega_1^2 r_1 + 2\zeta_1 \omega_1 \dot{r}_1 + (A_{111}^{(1)} + h_{111}^{(1)}) r_1^3 + B_{111}^{(1)} r_1 \dot{r}_1^2 \\
+ (A_{212}^{(1)} + A_{122}^{(1)} + h_{122}^{(1)}) r_1 r_2^2 + B_{122}^{(1)} r_1 \dot{r}_2^2 \\
+ B_{212}^{(1)} r_2 \dot{r}_1 \dot{r}_2 + C_{111}^{(1)} r_1^2 \dot{r}_1 + (C_{122}^{(1)} + C_{212}^{(1)}) r_1 r_2 \dot{r}_2 \\
+ C_{221}^{(1)} r_2^2 \dot{r}_1 = \hat{f} \cos(\omega t),
\end{aligned}
\quad (12a) $$

$$ 
\begin{aligned}
\ddot{r}_2 + \omega_2^2 r_2 + 2\zeta_2 \omega_2 \dot{r}_2 + (A_{222}^{(2)} + h_{222}^{(2)}) r_2^3 + B_{222}^{2} r_2 \dot{r}_2^2 \\
+ (A_{112}^{(2)} + A_{211}^{(2)} + h_{112}^{(2)}) r_2 r_1^2 + B_{211}^{(2)} r_2 \dot{r}_1^2 \\
+ B_{112}^{(2)} r_1 \dot{r}_1 \dot{r}_2 + C_{222}^{(2)} r_2^2 \dot{r}_2 + (C_{121}^{(2)} + C_{211}^{(2)}) r_1 r_2 \dot{r}_1 \\
+ C_{112}^{(2)} r_1^2 \dot{r}_2 = 0,
\end{aligned}
\quad (12b) $$

where $h_{ipk}^{(j)}$ are the coefficients of cubic terms in equation (7), $A_{ipk}^{(j)}, B_{ipk}^{(j)}$ and $C_{ipk}^{(j)}$ arise from the cancellation of the quadratic terms, and are expressed by

$$ A_{ipk}^{(j)} = \sum_{l \ge i}^{S} g_{il}^{(j)} a_{pk}^{(l)} + \sum_{l=1}^{l \le i} g_{li}^{(j)} a_{pk}^{(l)}, \quad (13a) $$

$$ B_{ipk}^{(j)} = \sum_{l \ge i}^{S} g_{il}^{(j)} b_{pk}^{(l)} + \sum_{l=1}^{l \le i} g_{li}^{(j)} b_{pk}^{(l)}, \quad (13b) $$

$$ C_{ipk}^{(j)} = \sum_{l \ge i}^{S} g_{il}^{(j)} c_{pk}^{(l)} + \sum_{l=1}^{l \le i} g_{li}^{(j)} c_{pk}^{(l)}, \quad (13c) $$

where $g_{il}^{(j)}$ are the coefficients of quadratic terms in equation (7). Forcing term has been added at the end of the process and now appears in Eq. (12a). This is the second approximation used for building this ROM, as a time-invariant manifold is used. The most accurate solution would have consisted in constructing a periodically forced invariant manifold, see e.g. (Jiang et al. 2005b). However, this results in a very complicated formulation and time-consuming numerical calculations for constructing the ROM. The proposed method has the advantage of the simplicity, the quickness of the computation, and allows deriving a differential model that could be used easily for parametric studies. However, it is valid, strictly speaking, only for small values of $\hat{f}$.

With the NNMS method, the original 16 degree of freedom of the conventional Galerkin model have been reduced to two in equation (12). However, differently to POD method, the structure of the equations of motion is changed. In fact, quadratic nonlinear terms have been cancelled, but cubic terms involves both the transformed nondimensional displacement and the transformed nondimensional velocity.
## 3.4 DISCUSSION

After presentation of the two reduction methods, a first discussion on their theoretical settings is here provided in order to underline their abilities and limitations.

The POD method, which consists in finding the best orthogonal hyper-planes that contains most information, is essentially a linear method. This can be seen as an advantage since few manipulations, involving linear algebra only, are needed to construct the ROM. The key formula of the method, Eq. (14), is an eigenvalue problem. On the other hand, the linear essence of the method may be a drawback, as curved subspaces are generally more suitable to capture clouds of points with complicated shapes. A NNM, being an invariant manifold in state space, is a curved subspace, so that the NNMs reduction method is essentially nonlinear. The invariance property is the key that allows finding the lowest dimensional subspaces that contains dynamical properties, since dynamical motions do not stay within any other subspace that do not share this invariance property. This is the main advantage of the NNMs method as compared to the POD. It is expected that fewer NNM are necessary than POD modes. This will be illustrated in sections.

The POD method is global in the sense that it is able to capture any motion in state space and furnishes the adapted basis for decomposing it. This is an advantage as compared to the asymptotic NNMs method used here, which relies on a local theory. The third-order development, Eq. (11), is valid for small values of the modal amplitudes. Moreover, the use of a time-invariant manifold prevents for reliable results for large values of the amplitude of external forcing. When increasing the non-linearities by feeding more energy into the system, the results provided by the asymptotic NNMs method are awaited to deteriorate. This is not the case for the POD, if one has taken care to construct its POD-based ROM with clouds of points that are significant for a large range of values of the nonlinearity. In this context, it has already been argued (Kerschen 2003, Amabili et al. 2006) that a chaotic response is the best candidate for building the POD.

Finally, the two methods differ radically in the way the ROM is built. For the POD, one must mandatory have a response of the system to be able to build the ROM. In the present context of a completely theoretical model, this is a drawback since one must compute time responses to be in position of reducing. Moreover, it has been underlined by Amabili et al. (2003, 2006) that the choice of these time responses is not an easy task that could not be done blindly. By comparison, the asymptotic NNMs method do not need any response of the system, but dynamical properties only, that are here provided by Eq. (7), i.e. after projection of the PDE with the Galerkin method. With at hand the eigenvalues of the linear part (eigenfrequencies $\omega_j$ and damping coefficients $\zeta_j$) and the nonlinear coefficients $g_{il}^{(j)}$ and $h_{pk}^{(j)}$, the nonlinear change of coordinates, Eq. (11), can be applied directly to obtain the ROM. As the coefficients in Eq. (11) are computed once and

for all, application of the method is easy and do not ask for computation time. In the next section, all these conclusions will be illustrated with the numerical results.

## 4. NUMERICAL RESULTS

The simply supported, water-filled circular cylindrical shell (without imperfections) investigated by Amabili (2003) is considered, with the following dimensions and material properties: $L = 520$ mm, $R = 149.4$ mm, $h = 0.519$ mm, $E = 2.06 \times 10^{11}$ Pa, $\rho = 7800$ kg/m$^3$, $\rho_F = 1000$ kg/m$^3$ and $\nu = 0.3$.

<figure><img src="image_2.png" /><figcaption>FIGURE 1. Maximum amplitude of vibration versus excitation frequency, excitation 3 N; conventional Galerkin model, 16 dofs. (a) Maximum amplitude of $A_{1,n}(t)$, driven mode; (b) maximum amplitude of $B_{1,n}(t)$, companion mode; 1, branch “1”; 2, branch “2”; BP, pitchfork bifurcation; TR, Neimark-Sacker (torus) bifurcations. —, stable periodic solutions; ■, stable quasi-periodic solutions; —, unstable periodic solutions.</figcaption></figure>
Numerical calculations have been performed for the fundamental mode ($m = 1, n = 5$) of the water-filled shell. The natural frequency $\omega_{1,n}$ of this mode is 79.21 Hz, according to Donnell's theory of shells; modal damping $\zeta_{1,n} = 0.0017$ is assumed. The two reduction methods will be compared on their ability to recover frequency-response curves and time responses, for moderate values of the amplitude of forcing, respectively 3 N and 8 N.

The response of the fundamental mode of the water-filled shell to harmonic point excitation of 3 N at $\tilde{x} = L/2$ and $\tilde{\theta} = 0$ has been computed by using the conventional Galerkin model with 16 degrees of freedom (dofs); result is given in Figure 1. The solution presents a main branch “1” corresponding to zero amplitude of the companion mode $B_{1,n}(t)$; this branch has pitchfork bifurcations (BP) at $\omega/\omega_{1,n} = 0.9714$ and at 1.0018, where branch “2” appears. This new branch corresponds to participation of both $A_{1,n}(t)$ and $B_{1,n}(t)$, giving a travelling wave response. Branch “2” undergoes two Neimark-Sacker (torus) bifurcations (TR), at $\omega/\omega_{1,n} = 0.9716$ and 0.9949. Amplitude-modulated (quasi-periodic) response is indicated in Figure 1 for $0.9716 < \omega/\omega_{1,n} < 0.9949$, i.e. bracketed by the two Neimark-Sacker bifurcations.

<figure><img src="image_3.png" /><figcaption>FIGURE 2. Time response at excitation frequency $\omega/\omega_{1,n} = 0.991$, excitation 3 N; conventional Galerkin model, 16 dofs: (a) modal coordinate $A_{1,n}(t)$ associated to the driven mode; (b) modal coordinate $B_{1,n}(t)$ associated to the companion mode.</figcaption></figure>

The quasi-periodic time response of the shell for excitation of 3 N at frequency $\omega = 0.991\omega_{1,n}$, branch “2”, is reported in Figure 2 for the most significant generalized coordinates. This time response, which is more suitable than simple periodic responses to construct accurate POD models (Amabili et al.

2003), has been used to build a POD model. Both travelling waves around the shell in opposite directions have been taken into account (one obtained by direct integration and the other by changing the sign to the generalized coordinates associated to $\sin(n\theta)$ terms in equation (8)) to construct the POD reduced-order model. The optimal number of proper orthogonal modes $\tilde{K}$ to be retained in the reduced-order model can be estimated by plotting $\sum_{i=1}^{K} \lambda_i / \sum_{i=1}^{K} \lambda_i$ as a function of $\tilde{K}$; three proper orthogonal modes absorb practically all of the shell energy for the response at $\omega/\omega_{1,n} = 0.991$; therefore three proper orthogonal modes are used in the POD model. The first proper orthogonal mode is the driven mode, the second is the companion mode and the third is the axisymmetric mode.

Responses obtained by using the conventional Galerkin method (16 dofs) and by using the POD method (3 dofs) compare very well for excitation of 3 N, as shown in Figure 3 for both driven and companion modes; the main difference is a slight shift on the right of the first bifurcation point of branch “1”. It can also be observed that the natural frequency computed with the POD model is practically identical to the one computed with the Galerkin model. Figure 3 also shows the response computed with the NNMs method with only two degrees of freedom. It can be observed that also the response computed with the NNMs method compare very well with the original Galerkin model, with the curves just very slightly shifted in the left and with exact qualitative behaviour. In this case the maximum vibration amplitudes reach about 1.5 h for the driven mode and 0.9 h for the companion mode.

In order to compare results also in time domain, the quasi-periodic responses ($\omega/\omega_{1,n} = 0.991$, excitation 3 N) computed with the POD and NNMs models are reported in Figures 4 and 5, respectively; this response is more critical to be reproduced by reduced-order models than simple periodic responses. Whereas the response computed with the POD is in reasonably good agreement with the one in Figure 2 obtained with the Galerkin model, the response calculated by using the NNMs model is practically coincident with this.

Figure 6 has been obtained with the same three models (Galerkin, NNMs and POD for which the same equations obtained with time response for excitation of 3 N has been used), but for excitation of 8 N. In this case the maximum vibration amplitudes reach about 3 h for the driven mode and 2.5 h for the companion mode; differences among the three models become much more significant than in the previous case. In particular, the POD model is relatively close to the original Galerkin model, the main difference being the first bifurcation point of branch “1”, which is now significantly shifted on the right giving rise to a significant difference in the qualitative behaviour of the two models. The NNMs model has qualitatively the same behaviour of the original Galerkin model, but the response is significantly shifted on the left giving rise to model overestimating the softening nonlinearity of the system.
<figure><img src="image_4.png" /></figure>

<figure><img src="image_5.png" /><figcaption>FIGURE 3. Maximum amplitude of vibration versus excitation frequency, excitation 3 N; conventional Galerkin model, POD model with 3 modes and NNMs model with 2 modes. (a) Maximum amplitude of $A_{1,n}(t)$, driven mode; (b) maximum amplitude of $B_{1,n}(t)$, companion mode. ——, conventional Galerkin model (16 dofs); _____, POD model (3 dofs); - - -, NNMs model (2 dofs).</figcaption></figure>

It can be observed here that the POD model could be improved by using a time response computed for excitation of 8 N to find the proper orthogonal modes; however it is interesting here to investigate the robustness of a reduced order model to changes in the system parameters, therefore it is convenient to use the same model. On the other hand, the NNMs model is built once and for all, and may not been changed when varying the amplitude of the forcing. The observed differences with the reference solution are the consequences of the two approximations used to build it: asymptotic development and time-invariant manifold.

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

<figure><img src="image_7.png" /><figcaption>FIGURE 4. Time response at excitation frequency $\omega / \omega_{1,n} = 0.991$, excitation 3 N; POD model, 3 dofs: (a) modal coordinate $A_{1,n}(t)$ associated to the driven mode; (b) modal coordinate $B_{1,n}(t)$ associated to the companion mode.</figcaption></figure>

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

<figure><img src="image_9.png" /><figcaption>FIGURE 5. Time response at excitation frequency $\omega / \omega_{1,n} = 0.991$, excitation 3 N; NNMs model, 2 dofs: (a) modal coordinate $A_{1,n}(t)$ associated to the driven mode; (b) modal coordinate $B_{1,n}(t)$ associated to the companion mode.</figcaption></figure>
<figure><img src="image_10.png" /><figcaption>FIGURE 6. Maximum amplitude of vibration versus excitation frequency, excitation 8 N; conventional Galerkin model, POD model with 3 modes and NNMs model with 2 modes. (a) Maximum amplitude of $A_{1,n}(t)$, driven mode; (b) maximum amplitude of $B_{1,n}(t)$, companion mode. ——, conventional Galerkin model (16 dofs); _____, POD model (3 dofs); - -, NNMs model (2 dofs).</figcaption></figure>

## 5. CONCLUSIONS

Results show that both the proper orthogonal decomposition (POD) and the nonlinear normal modes (NNMs) methods are suitable to build reduced order model of a water-filled shell. In particular, a larger reduction of the model is possible by using the NNMs method. However, the asymptotic formulation used here in the NNMs method does not make it suitable to study very large vibration amplitudes, where the POD model performs better.

These results have to be related to the properties of the two methods. The nonlinear character of the invariant manifolds

that defines NNMs allows better reduction than the POD method, which is a linear decomposition. On the other hand, the global nature of the POD provides very robust results even for complex dynamics. However, it has been found that, although being an approximation, the qualitative behaviour of the NNMs compares very well with the original solution, until the validity limits are attained. It has been found numerically that these limits are not that small: amplitude of vibration up to $3h$.

Construction of the NNM-based reduced-order model with the asymptotic method is direct and does not need for intensive computations, as a single nonlinear change of coordinates, computed once and for all, is required. The method can be thus blindly applied, provided the Galerkin projection has been performed on a large number of modes. On the other hand, for the POD, particular care must be placed in the choice of the time responses used to build it, as already discussed in Amabili et al. (2003).

To conclude, the investigations conducted here show that for moderate vibration amplitude, the asymptotic NNMs method provides more reduced equations that always recover the qualitative behaviour. However, the method must be modified in order to bypass its main limitation, which is due to the asymptotic development. Unfortunately, only numerical solutions are possible, thus leading to an intense numerical effort in order to build the reduced-order model. Hence, for very large vibration amplitude and large range of parameter variations, the POD method still performs better due to its global nature.

## ACKNOWLEDGMENTS

The first author (M.A.) gratefully acknowledges fundamental cooperation on related topics with Francesco Pellicano, Michael P. Païdoussis and Abhijit Sarkar. The second author (C.T.), gratefully acknowledges Gérard Iooss for stimulating discussions on the mathematical aspects of the problems discussed herein.

## REFERENCES

* Amabili, M., 2003. Theory and experiments for large-amplitude vibrations of empty and fluid-filled circular cylindrical shells with imperfections. *Journal of Sound and Vibration* 262, 921-975.

* Amabili, M., Païdoussis, M. P., 2003. Review of studies on geometrically nonlinear vibrations and dynamics of circular cylindrical shells and panels, with and without fluid-structure interaction. *Applied Mechanics Reviews* 56(4), 349-381.

* Amabili, M., Sarkar, A., Païdoussis, M. P., 2003. Reduced-order models for nonlinear vibrations of cylindrical shells via the proper orthogonal decomposition method. *Journal of Fluids and Structures* 18, 227-250.

* Amabili, M., Sarkar, A. Païdoussis, M. P., 2006. Chaotic vibrations of circular cylindrical shells: Galerkin versus reduced-order models via the proper orthogonal decomposition method. *Journal of Sound and Vibration* 290, 736-762.

* Aubry, N. Holmes, P., Lumley, J. L., Stone, E., 1988. The dynamics of coherent structures in the wall region of a turbulent boundary layer. *Journal of Fluid Mechanics* 192, 115-173.
Azeez, M. F., Vakakis, A. F., 2001. Proper orthogonal decomposition (POD) of a class of vibroimpact oscillations. *Journal of Sound and Vibration* 240, 859-889.

Bellizzi, S., Bouc, R., 2005. A new formulation for the existence and calculation of nonlinear normal modes. *Journal of Sound and Vibration* 287, 545-569.

Breuer, K. S., Sirovich, L., 1991. The use of the Karhunen-Loève procedure for the calculation of linear eigenfunctions. *Journal of Computational Physics* 96, 277-296.

Carr, J., 1981. Applications of centre manifold theory. Springer-Verlag, New-York.

Doedel, E. J., Champneys, A. R., Fairgrieve, T. F., Kuznetsov, Y. A., Sandstede, B., Wang, X., 1998. AUTO 97: Continuation and Bifurcation Software for Ordinary Differential Equations (with HomCont). Concordia University, Montreal, Canada.

Elphick, C., Tirapegui, E., Brachet, M., Coullet, P., Iooss, G., 1987. A simple global characterization for normal forms of singular vector fields. *Physica D* 29, 95-127.

Guckenheimer, J., Holmes, P., 1983. Non-linear oscillations, dynamical systems and bifurcations of vector field. Springer, New-York.

Iooss, G., Adelmeyer, M., 1998. Topics in bifurcation theory. World Scientific, New-York, second edition.

Jézéquel, L., Lamarque, C. H., 1991. Analysis of non-linear dynamical systems by the normal form theory. *Journal of Sound and Vibration* 149, 429-459.

Jiang, D., Pierre, C., Shaw, S., 2005a. The construction of non-linear normal modes for systems with internal resonance. *International Journal of Non-linear Mechanics* 40, 729-746.

Jiang, D., Pierre, C., Shaw, S., 2005b. Nonlinear normal modes for vibratory systems under harmonic excitation. *Journal of Sound and Vibration* 288, 791-812.

Kerschen, G., Feeny, B. F., Golinval, J.-C., 2003. On the exploitation of chaos to build reduced-order models. *Computer Methods in Applied Mechanics and Engineering* 192, 1785-1795.

Kerschen, G., Golinval, J.-C., Vakakis, A. F., Bergman, L. A., 2005. The method of proper orthogonal decomposition for dynamical characterization and order reduction of mechanical systems: an overview. *Nonlinear Dynamics* 41, 147-169.

King, M. E., Vakakis, A. F., 1994. Energy-based formulation for computing nonlinear normal modes in undamped continuous systems. *Journal of Vibration and Acoustics* 116, 332-340.

Lacarbonara, W., Rega, G., Nayfeh, A. H., 2003. Resonant non-linear normal modes. Part I: analytical treatment for structural one-dimensional systems. *International Journal of Nonlinear Mechanics* 38, 851-872.

Mikhlin, Yu. V., 1995. Matching of local expansions in the theory of non-linear vibrations. *Journal of Sound and Vibration* 182, 577-588.

Pellicano, F., Amabili, M., Païdoussis, M. P., 2002. Effect of the geometry on the non-linear vibration of circular cylindrical shells. *International Journal of Non-Linear Mechanics* 37, 1181-1198.

PesHECK, E., Pierre, C., Shaw, S., 2002. A new Galerkin-based approach for accurate non-linear normal modes through invariant manifolds. *Journal of Sound and Vibration* 249, 971-993.

Poincaré, H., 1892. Les méthodes nouvelles de la mécanique céleste. Gauthiers-Villars, Paris.

Rosenberg, R. M., 1966. On non-linear vibrations of systems with many degrees of freedom. *Advances in Applied Mechanics* 9, 155-242.

Sarkar, A., Païdoussis, M. P., 2003. A compact limit-cycle oscillation model of a cantilever conveying fluid. *Journal of Fluids and Structures* 17, 525-539.

Sarkar, A., Païdoussis, M. P., 2004. A cantilever conveying fluid: coherent modes versus beam modes. *International Journal of Non-Linear Mechanics* 39, 467-481.

Shaw, S., Pierre, C., 1991. Non-linear normal modes and invariant manifolds. *Journal of Sound and Vibration* 150, 170-173.

Shaw, S. W., Pierre, C., 1993. Normal modes for non-linear vibratory systems. *Journal of Sound and Vibration* 164, 85-124.

Sirovich, L., 1987. Turbulence and dynamics of coherent structures, Part I: coherent structures. *Quarterly of Applied Mathematics* 45, 561-571.

Slater, J. C., 1996. A numerical method for determining nonlinear normal modes. *Nonlinear Dynamics* 10, 19-30.

Touzé, C., Amabili, M., 2005. Non-linear normal modes for damped geometrically non-linear systems: application to reduced-order modeling of harmonically forced structures. *Journal of Sound and Vibration* (submitted).

Touzé, C., Thomas, O., 2006. Non-linear behaviour of free-edge shallow spherical shells: Effect of the geometry. *International Journal of Non-linear Mechanics* (accepted).

Touzé, C., Thomas, O., Chaigne, A., 2004. Hardening/softening behaviour in nonlinear oscillations of structural systems using non-linear normal modes. *Journal of Sound and Vibration* 273, 77-101.

Vakakis, A. F., Manevich, L. I., Mikhlin, Yu. V., Philipchuck, V. N., Zevin, A. A., 1996. Normal modes and localization in non-linear systems. Wiley, New-York.

Wolfram, S., 1999. *The Mathematica Book*, 4th edition. Cambridge University Press, Cambridge, UK.

Zahorian, S. A., Rothenberg, M., 1981. Principal component analysis for low-redundancy encoding of speech spectra. *Journal of the Acoustical Society of America* 69, 519-524.