Forced vibrations of circular plates: from periodic to chaotic motions

A numerical study of the transition from periodic to chaotic motions in forced vibrations of circular plates, is propose d. A pointwise harmonic forcing of constant excitation frequen cy Ω and increasing values of the amplitude is considered. Perfe ct and imperfect circular plates with a free edge are studied wi thin the von Ḱarmán assumptions for large displacements (geometric non-linearity). The transition scenario is observed for di fferent excitation frequencies in the range of the first eigenfreque ncies of the plate. For perfect plate with no specific internal resona nce relationships, a direct transition to chaos is at hand. For imp erfect plate tuned so as to fulfill specific internal resonance relat ions, a coupling between internally resonant modes is first observ ed. The chaotic regime shows an attractor of large dimension, an d thus is studied within the framework of wave turbulence. INTRODUCTION Thin structures such as plates and shells can easily experience chaotic vibrations when subjected to intensive lo adings. Predicting the appearance of these chaotic behaviour s s of particular importance for fluid-structure interaction p roblems, civil engineering applications or aeroelastic studies. De spit the large documentation nowadays available on chaotic vibr ations [1–3], most of the studies deal with few degrees of free dom systems and the application of the classical route to cha s ∗Address all correspondence to this author. (Feigenbaum scenario, Ruelle-Takens-Newhouse quasiperi odic route, Pomeau-Manneville intermittent behaviour, ...) to mechanical systems. Experimental observations made on numer ous plates, shells, cymbals and gongs [4–6] revealed a gener ic route from periodic to chaotic vibrations in harmonically f orced thin structures involving a quasiperiodic state where ener gy is exchanged between internally resonant modes, before the appe arance of Fourier broadband spectrum, which is characteristi c of the aperiodic motion. From the numerical viewpoint, the transition from periodic to chaotic vibrations in continuous systems has been consid ered by Awrejcewicz, Kryskoet al in a series of papers [7–9], with various numerical methods (finite-difference and Bubn ovGalerkin method) and for plates, cylindrical shells, panel s and sector-type spherical shells. Amabili studied the transit ion o chaotic vibrations for circular cylindrical shells and dou blycurved panels [3, 10, 11], but only in the vicinity of the fund amental frequency. Comparisons of experimental and theoret ical results are provided by Murphy et al for a plate subjected to large-amplitude acoustical excitation [12]. Nagai, Mur ayama et al provided very detailed analytical and experimental result s for a shallow cylindrical panel, with or without a concentra ted mass [13,14]. From the theoretical viewpoint, the long-time behaviour an d the existence of global attractors and inertial manifolds i n von Kármán dynamical equations have been studied by Chueshov and Lasiecka for various type of damping laws (viscous, stru c1 Copyright c © 2010 by ASME tural and thermoelastic damping on the structure, or linear and nonlinear boundary dissipation) [15–17]. On the other hand , for very lightly damped plates without strong finite-size effec ts, the framework of wave turbulence can be applied to study the energy repartition through lengthscales [18,19]. The case of plates have been treated recently in [20], and experimental invest igations were reported in [21,22]. This paper presents a numerical investigation of the transi tion from periodic to chaotic motions in harmonically force d, perfect and imperfect circular plates with a free edge. Impe rfect plates are selected so as to tune the first eigenfrequencies, hence enforcing them to fulfill particular resonance relationshi ps. For perfect plates in the vicinity of the seven first eigenfreque nci s, a direct transition is observed. For imperfect plates, inte rnal resonances are first excited before the chaotic motion. The critical values of the force amplitude that are needed to att ain he chaotic state is carefully studied. Numerical results show that the imperfection significantly lowers this critical value, which is fully consistent with the experimental observations. Fi nally, the convergence of the Galerkin truncation are systematica lly adressed, with respect to the critical force, the Lyapunov e xponents and the Fourier spectra of vibration, leading to a disc us ion of the results obtained in the chaotic state, with regard to t he occurrence of wave turbulence or low-dimensional chaos. MODEL EQUATIONS Perfect and imperfect circular plates The dynamic analog of the von Kármán equations for thin plates are used to model the large-amplitude vibrations. A c ircular plate is considered, fig. 1 shows the top view, the polar co o dinates and a cross-section defining the geometric imperfec tion. The starting point of the present study is the non-dimension al equations for imperfect plates presented in [23]. They read : ∆∆w+ ẅ= ε[L(w,F)+L(w0,F)]−2μẇ+ p, (1a) ∆∆F = − 1 2 [L(w,w)+2L(w,w0)], (1b) wherew is the transverse displacement, w0 the geometric imperfection,F the Airy stress function, μ an ad-hoc viscous damping, and p represents the external forcing. L is a bilinear operator defined by:


INTRODUCTION
Thin structures such as plates and shells can easily experience chaotic vibrations when subjected to intensive loadings. Predicting the appearance of these chaotic behaviours is of particular importance for fluid-structure interaction problems, civil engineering applications or aeroelastic studies. Despite the large documentation nowadays available on chaotic vibrations [1][2][3], most of the studies deal with few degrees of freedom systems and the application of the classical route to chaos * Address all correspondence to this author.
(Feigenbaum scenario, Ruelle-Takens-Newhouse quasiperiodic route, Pomeau-Manneville intermittent behaviour, ...) to mechanical systems. Experimental observations made on numerous plates, shells, cymbals and gongs [4][5][6] revealed a generic route from periodic to chaotic vibrations in harmonically forced thin structures involving a quasiperiodic state where energy is exchanged between internally resonant modes, before the appearance of Fourier broadband spectrum, which is characteristic of the aperiodic motion.
From the numerical viewpoint, the transition from periodic to chaotic vibrations in continuous systems has been considered by Awrejcewicz, Krysko et al in a series of papers [7][8][9], with various numerical methods (finite-difference and Bubnov-Galerkin method) and for plates, cylindrical shells, panels and sector-type spherical shells. Amabili studied the transition to chaotic vibrations for circular cylindrical shells and doublycurved panels [3,10,11], but only in the vicinity of the fundamental frequency. Comparisons of experimental and theoretical results are provided by Murphy et al for a plate subjected to large-amplitude acoustical excitation [12]. Nagai, Murayama et al provided very detailed analytical and experimental results for a shallow cylindrical panel, with or without a concentrated mass [13,14].
From the theoretical viewpoint, the long-time behaviour and the existence of global attractors and inertial manifolds in von Kármán dynamical equations have been studied by Chueshov and Lasiecka for various type of damping laws (viscous, struc-tural and thermoelastic damping on the structure, or linear and nonlinear boundary dissipation) [15][16][17]. On the other hand, for very lightly damped plates without strong finite-size effects, the framework of wave turbulence can be applied to study the energy repartition through lengthscales [18,19]. The case of plates have been treated recently in [20], and experimental investigations were reported in [21,22]. This paper presents a numerical investigation of the transition from periodic to chaotic motions in harmonically forced, perfect and imperfect circular plates with a free edge. Imperfect plates are selected so as to tune the first eigenfrequencies, hence enforcing them to fulfill particular resonance relationships. For perfect plates in the vicinity of the seven first eigenfrequencies, a direct transition is observed. For imperfect plates, internal resonances are first excited before the chaotic motion. The critical values of the force amplitude that are needed to attain the chaotic state is carefully studied. Numerical results show that the imperfection significantly lowers this critical value, which is fully consistent with the experimental observations. Finally, the convergence of the Galerkin truncation are systematically adressed, with respect to the critical force, the Lyapunov exponents and the Fourier spectra of vibration, leading to a discussion of the results obtained in the chaotic state, with regard to the occurrence of wave turbulence or low-dimensional chaos.

MODEL EQUATIONS Perfect and imperfect circular plates
The dynamic analog of the von Kármán equations for thin plates are used to model the large-amplitude vibrations. A circular plate is considered, fig. 1 shows the top view, the polar coordinates and a cross-section defining the geometric imperfection. The starting point of the present study is the non-dimensional equations for imperfect plates presented in [23]. They read: where w is the transverse displacement, w 0 the geometric imperfection, F the Airy stress function, µ an ad-hoc viscous damping, and p represents the external forcing. L is a bilinear operator defined by: When w 0 = 0 (no imperfection), one recovers the usual von Kármán equations for perfect plates [24]. In this contribution, all the quantities are non-dimensional, the relation with their dimensional counterpart are given in [23,24], where the interested reader can also find more details on the plate model and the freeedge boundary condition. The equations of motion (1) are discretized by expanding the unknowns w, w 0 and F onto the eigenmodes of the perfect plate, i.e. the modes of the linear system associated to Eqs. (1a,b) with w 0 = 0, via: where z g is the center of mass' offset, {a p } p=1...N 0 represent the projection of the imperfection onto the eigenmodes of the perfect plate, and q p , the modal amplitudes for the vibratory part, which will be our main unknown in the remainder of study. In the above equations, the {Φ i } i∈N are the transverse vibration mode shapes of the perfect plate and the {Ψ i } i∈N are membrane modes [23,24]. Assuming normalization of the modes and applying a usual Galerkin procedure, one finally obtains the following dynamical equations governing the evolution of the modal amplitudes, for all u ∈ [1, N]: Expressions for p u and Γ u rsp are: The cubic coefficients Γ u rsp are those of the perfect case, presented in [24]. The linear and quadratic coupling coefficients α u p and β u pr appearing in Eq. (4) stem from the geometric imperfection w 0 , and are thus expressed as functions of the amplitudes a p of the expansion of w 0 introduced in Eq. (3a). They write [23]:

Numerical details
In the remainder of the paper, a non-dimensional pointwise harmonic force of amplitude F and frequency Ω is considered, located at the edge of the plate so that (r 0 , θ 0 )=(1,0). The external distributed force reads: The angular location has been chosen so as to directly excite the cosine configuration, whereas the sine configuration of all the asymmetric modes is not linearly excited.
In order to simulate the response of lightly damped plates composed of metallic alloys, the damping coefficient has been set to µ = 0.002. All the quantities appearing in Eqs. (4) are analytic, computed and stored in subsequent arrays. The critical point in Galerkin expansions such as those assumed here in Eqs. (3) is the number of modes N and N F retained. N F controls accuracy of coefficents and is fixed at 12 [24]. Convergence with N is here studied by increasing N up to 35. Most of the results presented are computed with a 19-modes truncation, where asymmetric modes (2, 0) to (6, 0) are retained (resulting in 14 oscillator equations as each time sine and cosine configuration  Table 1. Non-dimensional frequencies of a perfect circular plate with a free edge, by increasing order. The first sixteen frequencies are given. To each asymmetric eigenfrequency (i.e. with k = 0) corresponds two companion modes, respectively sine and cosine configuration.
are present), and axisymmetric modes from (0, 1) to (0, 5). Table 1 recalls the first eigenfrequencies of a perfect plate. In the remainder, asymmetric companion modes will be distinguished with a letter C for cosine and S for sine mode, e.g. (2, 0,C). The resulting ordinary differential equations (ODEs) governing the dynamics are numerically integrated by using the DIV-PAG routine of the Fortran library IMSL, where a variable-order method based on Backward Differentiation Formulas (BDFs), also known as Gear's BDFs, is implemented. For analyzing the results when varying the amplitude of the external forcing, three usual techniques are used: stroboscopic Poincaré maps at the excitation frequency (with a zero phase) , power spectra of the transverse velocity and Lyapunov exponents. Figure 2 shows the result obtained by direct integration for an excitation frequency Ω=9.4 that is slightly over the first axisymmetric mode (0, 1), which is thus the directly excited mode in this simulation. Five other modal coordinates are shown, namely (2, 0,C), (2, 0, S), (3, 0,C) and (3, 0, S), and (1, 1,C). One can observe that in the periodic regime, q (0,1) has the most significant response. A branch switching is observed for F=0.9, corresponding to the loss of stability of the lower branch of the response in favour of the upper branch (a single Nonlinear Normal Mode (NNM) is excited thus a Duffing-type frequency response is expected for moderate values of the forcing). When increasing again the amplitude of the forcing, a sudden transition to a chaotic motion is observed for F=12. 19, where energy is present on all the generalized coordinates (only the first six are shown, but the model contains 19 modes, all of which being excited in the chaotic regime), and the amplitudes of vibration are very large.

PERFECT PLATE Generic result
The result shown in figure 2 is generic in the sense that numerous simulations have been performed for Ω varying between 1.5 and 25, and in quite all cases a direct transition from periodic to chaotic motion has been observed. Before the chaotic regime, non-zero (though small) values for modes (2,0,C) and (3,0,C) are observed. These are due to non-resonant couplings between the linear modes, that bends the non-linear normal mode (NNM) corresponding to the mode (0,1) that is directly excited. As the invariant manifold is curved and the periodic orbits computed are contained within it, small non-zero values are obtained for some modal amplitudes. However a strong coupling is not at hand, and the regime is periodic along the NNM corresponding to (0,1). The numerical experiment reported here can be seen as a research of the stability limits of the family of periodic orbits that compose the NNMs continuing the first eigenfrequencies, when the amplitude of the forcing is increased.

Convergence study
A first convergence study is presented here, by inspecting the number of modes needed to predict accurately the critical force F cr needed to obtain the chaotic regime. Table 2 sums up the obtained results, for three different frequencies Ω tested, and for an increasing number of modes N retained in the truncation for the transverse displacement, see Eq. (3b).
The results show that the model with 19 modes capture the good value of F cr in the three cases. It also highlights the fact that a very limited number of modes are needed to predict the correct F cr for low-frequency excitations, this number increasing when the excitation frequency increases. Only 3 modes are necessary for Ω=5.3, 7 modes for Ω=9.4, and 16 for for Ω=21. However, once the chaotic state established, all the modal coordinates have a significant amplitude, as already noted in Fig. 2, which means that the energy is spread all over the modes. Hence a major effect of increasing the number of modes is to lower the response amplitude of the slow-flow modes. Fig. 3 shows the amplitude of mode (0,1) for the same experiment with Ω= 9.4, and conducted with N = 7 and N = 27. N = 7 is the first converged value, which can be seen by the fact that the jump from the low to the high-amplitude branch occurs for the same value F=0.9, and the chaotic motion is obtained for the same F cr . The only difference resides in the amplitude of the response in the chaotic regime, which is larger for N=7 because all the energy can be spread between these 7 modes only. For N=27, a flux of energy allows excitation of all the retained modes; consequently the amplitude of the first are lowered. This result indicates that the numerical solution is probably not fully converged in the chaotic regime. This will be more deeply investigated in the next subsection, where the wave turbulence framework will allow for a better understanding of the chaotic regime. Figure 4 shows the complete bifurcation diagram for a perfect circular plate for excitation frequencies Ω ∈ [1. 5,25]. Each point in the figure corresponds to a numerical experiment as shown in fig. 2 for Ω =9.4. The generic observation is a direct transition from periodic to chaotic response, where the susbequent values of the critical force F cr needed to attain chaos are reported, which is logical as no internal resonance relationships are present between the first eigenfrequencies of the plate. Also note that as only a cubic non-linearity is at hand, only third-order internal resonance could destabilize the periodic responses, and they are more difficult to obtain. Figure 4 has been obtained with the 19-modes model until Ω=16, then a 27-modes model has been used. Values of the eigenfrequencies, as well as subharmonics and superhamornics of order three are also plotted.

Lyapunov and Power spectra
The chaotic regime is now more deeply investigated by first inspecting the Lyapunov exponents of the full system in the chaotic range. A simulation with Ω=9.4 is analyzed, with a 19modes model, giving birth to 38 Lyapunov exponents. Figure 5 shows the whole Lyapunov spectrum just after the onset of the chaotic regime, for F=13, computed as a by-product of the integration procedure, following the technique explained in [3,11]. The most striking result is the number of positive exponents: half of them are positive, half of them are negative. This indicates, by applying the Lyapunov dimension assumption, that the attractor have a dimension equal to that of the phase space. Hence a lowdimensional chaos is not at hand and the trajectories explore the whole phase space.
In order to get insight into the chaotic regime, the wave turbulence theory can be applied. Wave (or weak) turbulence (WT) provides a unified framework for out-of-equilibrium systems, provided the non-linearity is weak, and no intermittencies are at hand. The theory has been applied succesfully to plate dy-  namics with von Kármán assumptions to predict the shape of the power spectra [20]. Translated in the frequency domain, the prediction is that of a power spectra of the transverse velocity Pẇ( f ) which is independent of the frequency f . Different experiments led on very thin plates have revealed a discrepancy on this prediction, showing independently on two set-ups a clear dependence on the frequency for the velocity power spectra [21,22]. On the other hand, intermittency are not present and the persistence of waves is experimentally assessed, underlining that the main assumptions of the wave turbulence are fulfilled. The discrepancy is attributed to the presence of damping in the plate dynamics, which is not taken into account in the analytical results that assume a window of transparency and thus Hamiltonian dynamics. Figure 6 shows the power spectra of the transverse velocity obtained numerically with the model presented here; for a forc-  ing with Ω=5.3 (just below the first eigenfrequency) and F=20 (sufficient to have WT regime), and for three different truncations with an increasing number of modes, respectively N=10, 19 and 35. One can clearly see a first region, corresponding to the cascade regime where energy flux is present, where the theoretical prediction (a flat spectrum, independent of the frequency) is fulfilled. This region broadens from the first eigenfrequency where the injection of energy is set, to the higher mode that is contained in the truncation. Two peaks are then visible (at f = 25 and 39 Hz): they correspond to the two axisymmetric modes that are added to the truncation for stability reasons. This indicates that the more and more modes are added to the truncation, the more and more larger this cascade regime will be. Hence the model is not converged in the chaotic region. These results also highlights the fact that the convergence will be very difficult to obtain numerically : one has to include more and more modes until the dissipation scale. As the numerical burden of the Galerkin-type models used here increased very rapidly with the number of retained modes, these results show that another numerical method should be more efficient to simulate satisfactorily the WT regime.

IMPERFECT PLATE
An imperfect plate is now studied. The selected imperfection has the shape of the first axisymmetric mode (0,1), with a nondimensional amplitude of 0.45 (as the space variables are made nondimensional with respect to the thickness h, this means  Table 3. Non-dimensional eigenfrequencies of an imperfect circular plate, imperfection having the shape of mode (0,1) with an amplitude of an amplitude imperfection of half the thickness). The first eight eigenfrequencies of this imperfect plate are shown in Table 3. This amplitude of imperfection has been selected so as to fulfill the following internal resonance relationships: Hence the 2:1 internal resonance relationships are awaited to be easily activated when forcing around the involved eigenfrequencies, thus leading more easily to a coupled regime involving more than one NNM. This is confirmed by the numerical simulations. Two results are shown for Ω=10.2 and Ω=10.6, respectively in Fig. 7 and 8, with a model composed of 17 modes. The excitation frequency has been selected around the resonance frequency of mode (0,1), thus mode couping with (2,0) and/or (1,1) is awaited.
For Ω=10.2, one can observe that the 2:1 internal resonance with mode (2,0) only is activated, for a very small value of the forcing amplitude, namely F=0.81. As already pointed out for shells [25] for this case of 1:1:2 internal resonance (the 1:1 resonance arising from the existence of the two companion modes, here namely (2,0,C) and (2,0,S)), the coupling exist with only one companion mode, simultaneous coupling of the two being impossible. This is confirmed here on the transition scenario where the exchange of energy is observed to mode (2,0,S) only. The coupling is also characterized by a period-doubling in the temporal response which is clearly visible on the Poincaré sections with the appearance of two points. On the other hand, the other 2:1 relationship with mode (1,1) is not activated. Finally the chaotic regime is obtained for a small value of the forcing, F=2.81. These numerical results underlines the fact that adding imperfection leads to create more easily internal resonance relationships through quadratic nonlinearity, and is clearly a way to enhance the occurrence of chaotic vibrations for small values of external forcing. Figure 8 shows the result for a slightly higher value of the excitation frequency, Ω = 10.6. In this case the scenario before the chaotic regime is little bit more complicated. In the first region for F ∈ [0, 0.8], the coupling with mode (2,0) is immediately activated and is directed to mode (2,0,C). Residual nonzero value on mode (2,0,S) with decreasing amplitude is also observed which is due to a transitory response still being present in the Poincaré section. Amplitude modulation is then observed around F=0.8 ; this branch of solution being then lost in favour of a clear coupling with mode (2,0,S) for F ∈ [0.83, 2.2]. At F=2.2, the other 2:1 relationship is activated and mode (1,1,C) comes into play with a non-zero amplitude. This bifurcation leads to a coupling between mode (0,1) and (2,0) changing to the cosine mode (2,0,C). Finally, the wave turbulence regime is obtained for F=6.33.
As a conclusion for the imperfect plate, numerical simulations evidenced the fact that the imperfection favours the apparition of the chaotic regime for smaller values of the forcing amplitude, a fact that is mainly due to the appearance of quadratic non-linearity which, in turn, makes possible second-order internal resonance relationships. When these resonance relationships are fulfilled, a first instability occurs in the system where the directly forced NNM loses stability in favour of a coupled regime with all the modes involved in internal resonance. Then, for significantly smaller values, the wave turbulence regime is obtained. Other imperfections are considered in [26] as well as a brief review of experimental results, which definitely assess the transi-   tion scenario from periodic to chaotic motion.
Finally Fig. 9 shows the power spectra of the transverse velocity for the perfect plate and two imperfect plates with an increasing amplitude, respectively 0.45 and 1.86 times the thickness. An open question in the wave turbulence regime is the effect of the imperfection on the predicted power spectra, the theoretical prediction [20] being for perfect plates only. The presence of quadratic nonlinearity enforces a 3-waves interaction process in the cascade regime, which could significantly change the frequency dependence in the spectra. Moreover, experimental measurements shows that 3-waves interactions seem to be predominant through the power-law dependence of the spectra with respect to the injected power [27]. The spectra reported in Fig. 9 have been obtained with a 30-modes model, they clearly evidenced that the power-law dependence of the spectra appears to be unaffected by the presence of an imperfection. Further research in this topic will be to consider realistic damping laws, as well as to compute, from the numerics, the dependence with the injected power, in order to get a full comparison between theory and experiments.

CONCLUSION
The transition from periodic to chaotic motion in forced vibrations of plates has been numerically studied. The whole scenario can be interpreted as a transition to turbulence in a solid system. It involves either a direct transition, if no internal resonance relationships exist between the eigenmodes of the plates, or a generic quasiperiodic motion, when internal resonance are present, where energy is first spread to those internally resonant modes. This quasiperiodic state can degenerate to a periodic regime when the resonance relationships are the simplest ones, i.e. of the form 2ω p ≃ ω j (quadratic nonlinearity), or 3ω p ≃ ω j (cubic nonlinearity). The chaotic regime can be interpreted in the framework of wave turbulence. Numerical simulations allows recovering the analytical prediction. Discrepancies with the measured spectra revealed in [27] asks however for further refined numerical simulations, taking into account a realistic damping.