Multistability and Modal Interactions in Periodic 2D Coupled Pendulums Array

The collective dynamics of an array of periodic two dimensional (2D) coupled pendulums under harmonic horizontal base excitation is investigated. The coupled differential equations governing the nonlinear vibrations of the considered system have been solved using an analytical-numerical solving procedure, based on the multiple scales method coupled with standing wave decomposition. It allows the identiﬁcation of complex and wide variety of nonlinear phenomenon exhibited by the periodic nonlinear structure. The frequency responses for several coupled pendulums were calculated in order to analyze the stability, the modal interactions and the bifurcation topologies resulting from the collective dynamics of the coupled pendulums, while high-lighting the large number of multimodal solutions for a small number of coupled pendulums. The complexity and the multival-udness of the responses were illustrated by a study of basins of attraction which display the large distribution of the multi-mode branches.


INTRODUCTION
Recent researches focused on the study of wave propagation in nonlinear periodic media seen the existence of complex wave phenomena such as multi-modal solutions, modal interactions, solitary waves and discrete breathers in addition to those exhibited by linear periodic structures.A large number of phys-ical systems in different scientific areas such as materials, optics, chemistry and acoustics can be modeled as one-dimensional (1D) periodic nonlinear chain.In addition, the localization phenomenon generated by the interaction between the nonlinearities and a discrete system has attracted a lot of attention in recent years, because of its important role in the qualification as well as quantification of systems operations.One of the most popular localization phenomena, that have attracted the interest of physicists, is the nonlinear energy localization, which are called intrinsic localized modes ILMs.For instance, Sivers et al. [1] proved theoretically that ILMs could occur in pure anharmonic crystals.Recently, in the field of engineering, such as microelectromechanical systems (MEMS), ILMs were observed when micromechanical oscillators arrays were excited at high frequencies [2], and their characteristics were theoretically investigated [3].
In mechanics, Several researches focused on ILMs in pendulums array, we can mention Ikeda et al. [4] that investigated the behavior of ILMs for an array of coupled pendulums subjected to horizontal or vertical excitations.In addition, modulation instability of the uniform mode in a pendulum array with harmonic coupling and horizontal, sinusoidal driving which at longer times makes possible the formation of stable ILMs has been observed [5].Moreover, an experiment using a coupledpendulum chain has been reported to demonstrate various aspects of wave motion [6].Recently, The influence of adding external harmonic excitation on the intrinsic localized modes of coupled 1 pendulums chains has been investigated [7].
Taking into account the considerable work in the literature devoted to analyze the richness of the nonlinear dynamical effects found in an array of coupled pendulums [8], 2D systems are expected to present additional dynamics phenomena not realizable in the 1D case.Recently, the first analysis of a 2D array of coupled pendulums proposed a computational model to study the nonlinear dynamics of such system under parametric excitation [9].Although the dynamics of coupled nonlinear pendulums was thoroughly investigated in the frequency and time-space domains, there is a real need to perform profound analysis of the collective dynamics of such systems in order to identify practical relations with the nonlinear energy localization phenomenon in terms of modal interactions and bifurcation topologies [10].
In this work, we investigate the collective dynamics of periodic 2D coupled pendulums array under external excitation.The coupled equations governing the nonlinear vibrations of the considered system have been solved using an analytical-numerical solving procedure, which allows a detailed study of complex model interactions of the periodic structure.The multiple scales method coupled with standing wave decomposition were used to transform the nonlinear system into a set of complex algebraic equations.The responses frequency for several coupled pendulums were plotted in the frequency domain in order to analyze the stability, the modal interactions and the bifurcation topologies resulting from the collective dynamics of the coupled pendulums, while highlighting the large number of multimodal solutions for a small number of coupled pendulums.The complexity and the multivaludness of the responses were illustrated by a detailed study of basins of attraction which display the large distribution of the multi-mode branches.

THEORETICAL ANALYSIS Problem Formulation
Figure 1 shows the discrete model of 2D coupled N × M -pendulums subjected to harmonic horizontal base excitation.The system is composed of N equidistant axles A n .Each axle has M equidistant pendulum "i, j" attached perpendicularly to the axle and designed as a uniform rod of length l and mass m fixed at the end.φ i, j is its rotational degree of freedom (dof) from its stable equilibrium position.Pendula "i, j"and "i + 1, j", "i, j" and "i, j + 1" are connected respectively by torsional and linear springs k x and k y .The periodic system is subjected to harmonic horizontal base excitation y = y 0 cos(Ωτ), with fixed boundary conditions φ i,M+1 = φ i,0 = φ N+1, j = φ 0, j = 0.
The kinetic energy T of the system can be expressed as  Where (y i, j , z i, j ) represents the Cartesian coordinates of the pendulum {i, j}, given as and the dot denotes typically the derivative with respect to time.Substituting Eq. (2) into Eq.( 1) the kinetic energy becomes The potential energy of such a problem can be written as: where g represents the acceleration of gravity and ∆ j, j+1 the distance between pendulums {i, j} and {i, j + 1}, which can be written as ( The Lagrange equations of motion can be written in the following form with L the Lagrangian operator defined by L = T −V , Q i, j = αl φi, j represents the generalized forces applied to the pendulum {i, j} and α the viscous damping coefficient. Substituting Equations ( 3) and ( 4) into the Lagrange Equation (6), we obtain the equation of motion of each pendulum {i, j} By expanding sin(φ i, j ) and cos(φ i, j ) into Taylor series up to the third and first order respectively, the equations ( 7) can be written as: In the following sections, and since we are interested in localized modes phenomenon, the third order expansion related to the linear coupling will be neglected.

Dimensionless equations
For convenience and equation simplicity, the following nondimensional variables are introduced: where φ D = my 0 Ω 2 αω 0 is the dynamic displacement of the associated linear system, Q = ω 0 ml α and ω 0 = g l represents the natural frequency of pendulum "i, j".After dropping Eqs.(9) into Eq.( 8) and dividing by α m 2 l 2 y 0 ω 0 Ω 2 , we obtain the following nondimensional system of equations The aim of this paper is to present the influence of the modes localization on the collective dynamics of a 2D pendulums array.Therefore, a study on the natural frequency is required.

Linear study
Equation (10) describes a nonlinear Duffing system coupled by linear springs k x and k y in the x and y directions respectively and subjected to harmonic horizontal base excitation.To determine the natural frequencies and their eigenvectors Equation (10) can be written in matrix form as with . .] the N×M-dimensional variable vector, the excitation vector ) , F NL (θ ) is the nonlinear stiffness vector and , where It is convenient to express dimensionless eigenvalues ω I, j in the following form.For example, It is Convenient to express all natural frequencies relative to the same dimensional reference frequency which is 1, so that ω i, j = (1 + λ i ∆ x + λ j ∆ y ) (i = 1, . . ., N and j = 1, . . ., M), (13) where λ i and λ j are the coefficient depending on the mode i, j of vibration, ∆ x = k x mgl and ∆ y = k y l mg .We suppose that k x << mgl and k y << mg l .Consequently, ∆ x and ∆ y << 1 and: This assumption leads to the creation of linear closed modes which permits to study the effect of the mode localization on the collective dynamics.

Semi-analytical approach
Several methods have been developed to solve similar nonlinear differential systems (10).However, these methods are time consuming and not-suitable to solve a high number of nonlinear differential equations.In this paper, we use the multiple scales method [11] as an approximate analytical solution to solve nonlinear differential systems.The main advantages of the present approach are the capacity of handling weakly nonlinear systems, which permits to visualize all physical responses branches and their properties in terms of modal interactions and bifurcation topology transfer.
We introduce the parameters listed in Equation ( 15) and set the external frequency an amount εω 0 Ω D away from the resonant frequency, whereby they contribute to the equations of amplitude.
We calculate the responses of the pendulums array to harmonic, horizontal excitation using secular perturbation theory, by expanding θ i, j (t) of a sum of standing-wave modes with slowly varying amplitude [12].
T ) sin(nq r ) sin(md p )e iω r,p t + c.c.) where T = εt is a slow time variable, that authorizes the complex amplitude A r,p (T ) to vary slowly in time.Since we proposed fixed boundary conditions, the possible wave components q r , q p can be given as After replacing the proposed solution ( 16) into the equation of motion, we can get at the order of ε 3 2 , N × M equations of the form: i, j = ∑ r ∑ p (r, p th secular term)e iω r,p t + other terms We must eliminate the N × M secular terms so that θ i, j (1) remain finite which allow us to determine the equations for the slowly varying amplitudes A r,p (T ).To extract the equations for these amplitudes, we make use of the orthogonality of the modes, by multiplying the r, p th secular term by sin(nq r ) sin(md p ) and summing over n and m.We also express all normal frequencies relative to 1, so that: We find that the equation of the r, p th amplitude A r,p (T ) is given by: 2iω r,p dA r,p dT + iηω r,p A r,p + 2(γ x sin 2 ( q r 2 ) + γ y sin 2 ( In solving Equations (19), we write the periodic steady state solution of the form Substituting Equation ( 20) into (19), we obtain the required equation for the complex amplitudes a r,p .+2(γ x sin( q r 2 ) 2 + γ y sin(

TABLE 1. DESIGN PARAMETERS FOR THE CORRESPONDING PERIODIC STRUCTURE DEPICTED IN
This solving technique has been used by Lifshitz et al. [13] to calculate the response of nonlinear coupled microelectromechanical oscillators under parametric excitation and has been validated with respect to the experimental measurements of Buks and Roukes [14].The complexity of the nonlinear algebraic System (21), limits the possibility to solve it analytically for a high number of coupled pendulums.We may use the appropriate time integration procedure in order to solve the differential System (19) numerically, which implements a version of the fourth-order Runge Kutta integration method.

RESULTS AND DISCUSSIONS
In order to study the collective dynamics of the 2D periodic structure, two different configurations were proposed in Table 1, which highlights the influence of linear coupling on the non-linear behavior in terms of frequency range and vibrational amplitude.
For the case of 2 × 2 coupled pendulums (N = 2, M = 2) and for the first design parameters listed in Table 1, the natural dimensionless Eigen frequencies of the associated linear system are very close to 1.In addition, as we considered k y = k x l 2 , the linear springs are expected to act the same way.Figures 3 show the response intensities as a function of the detuning parameter Ω D .The Single Mode (SM) solution branches corresponds to the Trivial solution of the System (19) with A 1,2 = A 2,1 = A 2,2 = 0.The multimodal solutions DM and QM are generated by the coupling between two and all modes respectively.These curves were ploted to underline: first the effect of modes localization on the high number of solution branches for a given frequency and second the identical contribution of the non-excited modes on the responses.Although, the advantages of modes localization, considering small coupling parameters does not affect the systems rotational degree of freedom (dof).Figure 4, displays the rotational dof of 2 × 2 coupled pendulums under harmonic horizontal base excitation as a function of frequency, which corresponds to a single forced Duffing oscillator.In addition, all multimodal solution branches are distributed between resonant and non-resonant branch, hence its complicated to identify their contribution.
We decide to increase linear coupling parameters in such a manner to preserve the localization phenomena.Therefore, we consider the second design parameters listed in Table 1, so that the natural dimensionless Eigen frequencies of the associated linear system are ω 1,1 = 1.0331, ω 1,2 = 1.0336, ω 2,1 = 1.0989 and ω 2,2 = 1.0994.Figure 5 shows the response intensities of 2 × 2 coupled pendulums.Unlike to the first case, linear springs k x and k y play different roles.In particular, the multimodal solution branches resulting from the interaction between A 1,1 and A 2,1 or A 2,2 present larger frequency range and their bifurcation point B 2 exceeds the one that corresponds to the non-resonant branch of the SM (B 1 ).In addition, in regards to the systems rotational dof, multimodal solution branches present additional features distinct from those presented by a single forced Duffing oscillator (larger bandwidth and higher vibrational amplitude).
Taking into account the high number of solution branches founded above, for a given frequency, we choose to explore the practical appearance of each attractor according to their modal contribution type.Figure 7, displays the distribution of the basins of attraction of the multimode solution branches in the Nyquist plane (Re(a 11 ), Im(a 11 )) of the case of 2 × 3 coupled pendulums under harmonic horizontal base excitation, for the second design parameters listed in Table 1.For the configuration ω = 11.63, and random initial conditions Re(a 12 ) = 0.047, Re(a 13 ) = 0.2447, Re(a 21 ) = 0.1868, Re(a 22 ) = −0.1315,Re(a 23 ) = 0.2802, Im(a 12 ) = −0.2037,Im(a 13 ) = −0.3110,Im(a 21 ) = −0.3165,Im(a 22 ) = 0.1256 and Im(a 23 ) = −0.4189.The solutions DM, T M, QM, 5 th M and 6 th M are related to the branches generated by the interaction between two, three, four, five or all modes respectively.This distribution highlights the robustness of the multimodal solutions and the high contribution of the modes into the responses.

CONCLUSION
The collective nonlinear dynamics of a 2D periodic array of coupled pendulums was modeled under harmonic horizontal base excitation.An analytico-numerical method suitable for nonlinear localized modes problems has been considered, based on a perturbation technique, combined with standing wave decomposition, transforming the normalized system into a set of coupled complex differential equations.The effects of modes localization on the collective nonlinear dynamics have been studied for the case of 2 × 2 coupled pendulums.The resulting complex system has been numerically solved for two different configurations, highlighting the impact of linear coupling parameters on the complexity and the multivaludness of the responses.The robustness of these additional branches was illustrated by a study of the basins of attraction for 2 × 3 coupled pendulums.In practice, this model can be adjusted to increase the number of multimodal solutions, their configuration and their vibrational amplitudes.A complete quantitative study is a veritable challenge that we will attempt in a future work.

FIGURE 2 .
FIGURE 2. FRONT VIEW OF THE SYSTEM.

FIGURE 3 .
FIGURE 3. RESPONSE INTENSITY OF 2 × 2 COUPLED PENDULUMS AS A FUNCTION OF FREQUENCY, UNDER harmonic horizontal base excitation FOR THE FIRST SET OF DESIGN PARAMETERS LISTED IN TABLE 1.