N − 1 modal interactions of a three-degree-of-freedom system with cubic elastic nonlinearities

In this paper the N − 1 nonlinear modal interactions that occur in a nonlinear three-degree-of-freedom lumped mass system, where N = 3, are considered. The nonlinearity comes from springs with weakly nonlinear cubic terms. Here, the case where all the natural frequencies of the underlying linear system are close (i.e. ω n 1 : ω n 2 : ω n 3 ≈ 1 : 1 : 1) is considered. However, due to the symmetries of the system under consideration, only N − 1 modes interact. Depending on the sign and magnitude of the nonlinear stiffness parameters, the subsequent responses can be classiﬁed using backbone curves that represent the resonances of the underlying undamped, unforced system. These backbone curves, which we estimate ana-lytically, are then related to the forced response of the system around resonance in the frequency domain. The forced responses are computed using the continuation software AUTO-07p. A comparison of the results gives insights into the multi-modal interactions and shows the frequency response of the system to branches of the backbone curves that represent such interactions.


Introduction
Understanding the effects of modal interactions of coupled systems of nonlinear equations is an important step in being able to predict the subsequent dynamic response of the system. In this paper, we consider a three-degree-of-freedom (3-DoF) lumped mass system with cubic stiffness nonlinearities. In particular we consider the potential forced responses that can occur by analysing the backbone curves of the underlying undamped, unforced system. The justification of this approach lies in that the vast majority of engineering structures is characterised by low level of damping. This implies that their dynamical behaviour is largely determined by the dynamics of the associated Hamiltonian system [10].
The motivation for this study is the possibility for modes, in multi-degree-of-freedom nonlinear systems, to interact with each other [17]. These types of modal interaction have been previously studied because they are often related to unwanted vibration effects in engineering structures [16]. The majority of the existing literature is for undamped, unforced systems, and includes structural elements such as beams, cables, membranes, plates and shells-see for example [1,15,27]. Several different analytical and numerical approaches have been used to approach this type of problem, such as perturbation methods [18], nonlinear normal modes [13,22,23,25] or normal form analysis [2,11,21]. The majority of work in the literature on modal interaction is based on 2-DoF nonlinear systems where two modes interact, see for example [3,5,[8][9][10]13,29], although for continuous systems higher numbers of modes can typically be retained in the approximation, see for example [24]. Some work has been carried out on 3-DoF lumped mass systems in the context of nonlinear vibration suppression [12].
In a system with N modes, it is possible for response solutions to exist in which some ( i.e. N −i ≥ 2) and/or all N of the modes interact. In this paper, we consider the case where this type of N − i modal interaction can occur. More specifically, the case where N = 3 is analysed, i.e. only two ( i.e. N − 1) modes can interact. To do this we have chosen a specific configuration of an in-line 3-DoF nonlinear oscillator with small forcing and light damping. Due to the structural symmetry, one of the modes of this system is linear and the other two modes are nonlinear. So even though all the modal natural frequencies are close, the linear mode is not coupled with the other two. Hence the study of the 3-DoF system is reduced to the modal interaction of two coupled modes. Our analysis follows the approach developed by [3,9] to consider nonlinear modal interactions of a 2-DoF nonlinear oscillator. However, the effects of two modes interacting in this 3-DoF system has subtle differences from those described in the previous literature. This will be explained in Sect. 4.
The paper is structured as follows. In Sect. 2 we describe the configuration of this in-line oscillator. Then, we apply a normal form transformation method to the 3-DoF system to obtain its potential backbone curves expressions. In Sect. 3 the backbone curves of the system with the hardening nonlinearity are computed. These curves are then used to infer the dynamic behaviour of the system, which in turn can be used to interpret the forced, damped behaviour. In the Sect. 4, the backbone curves of the system with softening nonlinearity are presented. Lastly, the stability of the backbone branches are analysed and their relation with the forced response are shown. Conclusions are drawn in Sect. 5.

System description and analytical method application
The example system considered here is shown in Fig. 1. It consists of three lumped masses which are all grounded via linear springs, k, and viscous dampers, c. The mass in the middle links the two masses at two sides via the linear viscous dampers, c , and the cubic nonlinear springs. The elastic force characteristics of these nonlinear springs is F = k (Δx)+κ(Δx) 3 . Each of the masses is harmonically forced at the same frequency, Ω, and at the amplitude, P i . The equations of motion for this system are written in the form, where x,ẋ andẍ are the physical displacement, velocity and acceleration vectors, respectively. K is the linear stiffness matrix, P x is the external force amplitude vector and N x (x,ẋ) contains the nonlinear and damping terms which are assumed to be small. Their expressions are, and To proceed, the second-order normal form method is applied to approximately solve the equations of motion of this system. First, in order to decouple the linear terms, the linear modal transformation is applied to Eq. 1 to obtain the linear modal decomposition equations in terms of the new modal state q = {q 1 q 2 q 3 } T as, where Λ is a diagonal matrix of the squares of the corresponding linear natural frequencies ω 2 n1 = k/m, ω 2 n2 = (k + k )/m and ω 2 n3 = (k + 3k )/m. N q is the vector containing the damping and nonlinear coupling terms, where μ = κ/m, 2ω n1 ζ 1 = c/m, 2ω n2 ζ 2 = (c+c )/m and 2ω n3 ζ 3 = (c + 3c )/m. In the modal coordinates, the external force amplitude vector, P q = {P m1 P m2 P m3 } T , is where Φ is the matrix of the linear modeshapes. Here the modeshapes used for the linear transformation are {1, 1, 1} T , {1, 0, −1} T and {1, −2, 1} T for the three modes, respectively. Then, the forcing transformation, q → v, is supposed to be applied to remove the non-resonant forcing terms. Here we assume that the system is forced near resonance, this allows us to set q = v. Substituting this into Eq. where Lastly, applying the near-identity nonlinear transform v → u to Eq. 5 by using v = u + H(u,u) gives the equation, where P u = P q and N u = n u u * includes only the nonlinear resonant terms responding at the fundamental frequency of the corresponding modes. Here u and H(u,u) used in the transformation are the fundamental and harmonic components of v, respectively. Substituting v i → u i = u i p + u im into Eq. 3 in v gives the functional form of the u * vector and the coefficient matrix n v so that N v = n v u * as below, where the subscript i represents the i th linear mode. Then the homological equation is used to compute a matrix, β, for determining the resonance terms in the near-identity transformation [21]. The individual β i,l terms, can be computed using, where the subscript l denotes the l th term of u * , written as u * l and s npl and s nml represent the power indices of u np and u nm of the term u * l , respectively. Finally the variable, r , that describes the ratio between the response frequencies of mode 2 and 3, i.e. r = ω r 3 /ω r 2 , is introduced, so that for the example considered here β is, where η i = 2ω ni ζ i /u and in β, a dash indicates that the value is insignificant since the corresponding coefficient in n v is zero. The zero terms in β represent unconditionally resonant terms which have to be retained in N u . Furthermore, there are also additional conditionally resonant terms which are potentially set to zero depending on the value of r . For example, r = 1 will lead further zero terms in β. For all the terms for which β = 0, the corresponding terms in n v are set equal to those in n u . The modal interaction case considered here is when ω n1 : ω n2 : ω n3 ≈ 1 : 1 : 1. This allows us to assume all modal fundamental response frequencies are at the forcing frequency, i.e. Ω = ω r 1 = ω r 2 = ω r 3 , such that r = 1. Therefore, the terms in column 2, row 7 and 8 and column 3, row 13 and 14 in β become zero and the resulting dynamic equations of u i are, where U i and φ i are the fundamental response amplitude and phase of the i th mode, respectively, into Eq. 8 and balancing the coefficients of e jΩt and e − jΩt , we obtain the time-invariant equations for the forced response of the system, where p = e j2(|φ 2 −φ 3 |) . The |φ 2 − φ 3 | term represents the phase difference between mode 2 and 3.
To obtain the backbone curves, the unforced, undamped system needs to be considered. Therefore, by setting the damping and external force to be zero, i.e. ζ i = 0 and P mi = 0, in Eq. 9 gives Then successively setting (U 2 and U 3 ), (U 1 and U 3 ) and then (U 1 and U 2 ) to zero for Eq. 10 results in three single-mode backbone curve solutions labelled S1, S2 and S3 S1 : Further observing Eq. 10, it can be seen that mode 1 is linear and not coupled with the other two modes, while modes 2 and 3 could potentially interact with each other. This modal interaction is also affected by the value of p. If u 2 and u 3 are both present, i.e. U 2 = 0 and U 3 = 0, Eq. 10b and 10c can be rearranged and set equal, so that To ensure Eq. 14 is real, the phase difference terms need to be, p = ±1. Here, p = 1 and p = −1 represent the in-unison and out-of-unison resonances, respectively. A detailed discussion describing how the value chosen of p and the corresponding resonant types can be found in [3,9]. Firstly, by setting p = +1 yields two extra backbone curves, labelled S4 + and S4 − , with the phase difference (15) and their corresponding backbone branch expressions are same and given as, The case where p = −1 yields two further backbone curves, denoted S5 + and S5 − . They are characterised by the phase differences, Substituting p = −1 into Eq. 14 gives the response amplitude and frequency relationships, From inspection of Eqs. 16a, 18a, it can be seen that, since ω n3 > ω n2 , the sign of μ must be negative for the equations to have real solutions such that the backbone branches are physically realisable.

Backbone curves
When the nonlinear stiffness is positive, μ > 0, the solutions for Eqs. 16a and 18a are complex. Therefore, for the hardening case, the backbone branches S4 ± and S5 ± have no physical meaning and only S1, S2 and S3 exist. This means that there is no nonlinear modal interaction between the three backbone curves. Figure 2 shows the backbones curves for the system with the hardening nonlinearity where ω n1 = 1, ω n2 = 1.005, ω n3 = 1.0015 and the system parameters are m = 1, k = 1, k = 0.01 and κ = 0.05. All panels show the backbone curves in the projection of the response frequency against the displacement. The panels in the first row show the amplitude of the fundamental components of the modal response, u 1 , u 2 and u 3 , and the second row shows the amplitude of the fundamental response of the physical displacement of mass 1 and 2, x 1 and x 2 . Note that due to the system symmetry the backbone curve plot of mass 3 will be identical to that of mass 1 if it is shown. The S1, S2 and S3 branches are the single-mode backbone curves. It can be seen that S1 is linear and S2 and S3 show the characteristic shape of a hardening response.

Forced response
To show how the backbone curves can help facilitate the interpretation of the modal interaction of the nonlinear system, the forced response amplitude of three masses, X 1 , X 2 and X 3 , in the frequency domain for the hardening nonlinear system is shown in Fig. 3. To illustrate the relationship, the corresponding backbone curves in The forced response has been computed from an initial steady state solution, found with numerical integration in MATLAB, which is then continued in forcing frequency using the software AUTO-07p [4].
From Fig. 3, it can be seen that the response curve can be divided into three parts from left to right by two minima, -For the first resonance, the response curves of three masses are similar to that of the linear oscillator and are centred around S1. -For the second resonance, the familiar shape of the responses of a typical Duffing oscillator are following S2 and the jump phenomenon can be also observed on the right part of the curve for mass 1 and 3. Note that due to the special values of the linear modeshape of mode 2, i.e. {1, 0, −1} T , there is no peak for the response of mass 2 in Fig. 3b within this bandwidth. -For the third resonance, in Fig. 3a, the curve around S3 contains a loop where the upper trajectory is unstable which occurs from the addition and subtraction of the modal contributions. In Fig. 3b, c, the response of the typical Duffing oscillator can be observed enveloping S3.
For the first part of the response, mode 1 or linear mode, the resonant peak points are exactly on the linear backbone curve S1. For the other parts of the forced response curves, the curves cross the backbone curves S2 and S3 almost at the fold points. This is in line with the results found for systems with a lower number of degrees of freedom and shows how the backbone curves provide a good estimation of the position and

Backbone curves
For the softening case, μ < 0, Eqs. 16a and 18a have real solutions. Therefore, the in-unison, S4 ± , and outof-unison, S5 ± , resonant backbone curves are physi-cal. In Fig. 4, the backbone curves for the softening case where ω n1 = 1, ω n2 = 1.005, ω n3 = 1.015 and κ = −0.05 are shown. The first and second columns show the backbone curves of the modal states, u 1 , u 2 and u 3 , and the physical displacements, x 1 , x 2 and x 3 , respectively. As with the hardening case, the singlemode backbone curves, S1, S2 and S3, for mass 1 and where both mode 2 and 3 are activated, are of primary interest.
In Fig. 4, as expected the branch S1 is still a straight vertical line and S2 and S3 are curved, but with an opposite bending direction compared with those of the hardening system. Furthermore, it can be seen that S4 ± emanate from branch S2 and S5± bifurcate from branch S3. This type of bifurcation phenomenon has been previously studied in other systems, see for example [3,9]. A significant difference from the previous systems in the literature is that the S5 ± here is a vertical straight backbone curve which means that, like a linear mode, the resonant frequencies will not vary with the response amplitude. Also the backbone curve S5 ± has a finite length and it ends at S2, which implies that the out-of-unison resonance only happens within a certain amplitude range for this example.

Stability of the backbone curve
We now consider the stability of the backbone curves. Here, only the stability analysis of the backbone branch S2 is given in detail due to the fact that both S4 ± and S5 ± intersect with it, as shown in Fig. 4. Note that on backbone curve S2 which is the solution of Eq. 12, u 3 is equal to zero. So the stability of S2 can be determined by considering the dynamics of u 3 about its zero solution (Note that u 1 is not considered here due to its independence of u 2 and u 3 ). When the zero solution of u 3 is unstable, the S2 solution is also unstable. Now, the unforced, undamped equation of mode 3, Eq. 8c, is considered with setting c m3 = 0 and P m3 = 0, such that, As u 3 is considered around its zero solution, the u 3 p u 3m u 3 term is so small that it has been neglected. The stability of the system is found by considering the amplitude and phase of u 3 to be a slowly varying functions of time using the parameter ε to denote 'smallness'. A similar treatment can also be seen in [28]. Combined with u i = u i p + u im , we can write u 3 as Furthermore, the derivatives of u 3 can be written as, and where order ε 2 terms have been neglected. Substituting the expression ofü 3 into Eq. 19 and then balancing the coefficients of the e jω r 3 t and e − jω r 3 t terms gives, anḋ Now these equations can be expressed in the form, where U 3 = {U 3 p U 3m } T and the Jacobian matrix is, So the stability of the zero solution of u 3 can be assessed by considering the eigenvalues of the matrix, f U 3 , about the equilibrium solution U 3 = 0. So, the eigenvalues,  Note that due to the identical natural frequencies, the branches S5 ± overlap with S1, so the S5 ± backbone curves are indicated by the short cross lines for distinction λ, are given by, where U 2 2 = U 2 p U 2m has been used. From inspecting the above equation, it can be seen that the values of λ can be only be either purely imaginary or real. When all the eigenvalues, λ, are purely imaginary the system is marginally stable, while when λ are real the system is unstable. This implies that the bifurcation point on the backbone curve occurs when both of the roots of Eq. 25 are zero, such that, Then using the expression for S2 given in Eq. 12 gives the bifurcation points, and BP 2 : These are the same points where S4 ± and S5 ± intersect with S2, respectively. Note that since the linear natural frequencies of the modes are similar, the resonant response of the second and third modes will be at the same frequency. Hence, here, the relationship ω r 3 = ω r 2 = Ω has been used to simplify the equations.
Using the same approach, it can be also shown that the parts of S2 below bifurcation point BP 1 and above BP 2 are stable and the part between the two bifurcation points is unstable.
The same method was also used to predict that on S3 there is another bifurcation point, which is the intersection point of S3 and S5 ± . Also the stability condition of branch S3 is that the part below BP 3 is stable and the section above is unstable. All these bifurcation points are also shown in Fig. 4. Using bifurcation theory for Hamiltonian systems and comparing with other similar systems in the published literature, all the bifurcation points here are (Hamiltonian) pitchfork bifurcations [3,7,25].

Forced response
We now examine the relationship between the forced response and the backbone curves for the softening nonlinear case. From the backbone curve results, it can be seen that the modal interactions occur for this case. Therefore, a simple forced configuration where only mode 2 is excited is chosen as two kinds of interactions bifurcate from backbone curve S2. For this forcing configuration, three different force amplitudes are chosen which are [P m1 , P m2 , P m3 ] = [0, 0.4, 0]×10 −3 , [0, 0.9, 0]×10 −3 and [0, 1.5, 0]× 10 −3 and the damping ratio is chosen to be ζ = 0.001. Figure 5 shows the forced response results superimposed on the backbone curves for mass 1.
-For the small force amplitude situation, Fig. 5a, there is only one response curve which is centred around S2. This curve is the response of a typical softening Duffing oscillator and only composed the response of mode 2, u 2 . For this case, the force is insufficient to trigger the modal interaction or jump. -For the medium force amplitude situation, Fig. 5b, there are three response curves (1 blue and 2 green). The two green curves following S4 ± bifurcate from the single-mode response curve (blue one) at two secondary bifurcation points, respectively, and they are composed of the response of both mode 2, u 2 , and mode 3, u 3 . -For the large force amplitude situation, Fig. 5c, there are two additional response curves (the black ones) surrounding S5 ± which also bifurcates from the single-mode response curve of mode 2. On these two curves, both mode 2 and 3 are present as well.
From the results in Fig. 5, it can be noticed that for the situation where only one nonlinear mode is directly excited the magnitude of the force amplitude may affect the occurrence of the modal interaction. More specifically, some critical values of the force amplitude must be reached to result in the modal interaction. Furthermore, the mixed-mode responses following the inunison backbone curves would be triggered prior to those for the out-of-unison backbone curves when the external force amplitude increases.

Discussion of the force amplitude for triggering the modal interaction
Now, we consider at what forcing amplitude the modal interaction is triggered for the single-mode excitation situation. In Fig. 5b, c, it can be seen that the mixedmode response curves emanate from the single-mode response curve at the secondary bifurcation points. We infer from this that the modal interaction response curves will appear from the single-mode response curve when only one linear mode is forced [2,9]. So, the existence of modal interaction is studied by looking for bifurcation points on the single-mode response curve. The theory used for detecting the secondary bifurcation points is the same used for the stability of backbone curves in Sect. 4.2. When the zero solution for the response of the unforced mode is stable (unstable), the single-mode response curve of the other forced mode is also stable (unstable). As a result the point of neutral stability is the secondary bifurcation point. The case where only mode 2 is forced has been considered as an example here. Along the response curve composed of only u 2 (blue curves in Fig. 5), the modal coordinate u 3 is zero, so we are going to study the stability of the zero solution of u 3 [6]. Firstly, substituting trial solutions for u 3 ,u 3 andü 3 , Eqs. 20, 21 and 22, into Eq. 8c with P m3 = 0 (as only mode 2 is forced) and ω r 3 = Ω (as the fundamental components of all modal responses of this system are at the force frequency). Note that here an accuracy of order ε 1 is considered. As the damping is assumed to be small which allows us to regard ζ i to be order ε 1 , so the second expression part of Eq. 21 can be ignored after implying the damping coefficient. Then balancing e jΩt and e − jΩt , the equations of U 3 = {U 3 p U 3m } T can be obtained as, where Then the eigenvalues , λ, for F U 3 can be obtained from, wherê Therefore, the stability of the zero solution, U 3 = 0, can be determined from the eigenvalues, λ i , i = 1, 2.
The roots of Eq. 31 are given by, For all parameter values in this exampleâ andb are positive. So when one of the eigenvalues is zero and another is negative real, the system is neutrally stable. Therefore, the bifurcation occurs whenĉ is zero, such that, Combining Eq. 33 with the equation for the singlemode response curve for u 2 obtained by setting U 3 = 0 in Eq. 8(b) expressed as below, it allows us to find the position and number of their intersection points. To ensure the solutions are physically reasonable, only the intersection points at U > 0 and Ω > 0 are considered. Those points are secondary bifurcation points that can be used to predict the onset of the modal interaction. There are three possible cases: (1) if there is zero or one intersection point, there will be no modal interaction. Fig. 6a. (2) if there are two or three intersection points, the modal interaction response following the in-unison backbone curves will exist. Fig. 6b.  To make a comparison, the corresponding response curves (black) in Fig. 5 are also shown here. From Fig. 6, it can be seen that the singlemode response curve results approximated using the normal form technique and those calculated by AUTO are in good agreement. In addition, the intersections of Eqs. 33 and 34 are very close to the secondary bifurcation points for every case. The results show that the analytical method presented in this paper can be applied to detect the occurrence of modal interaction and identify the position of bifurcation points. Furthermore, when only mode 3 is excited, the same analysis method can be applied to find the bifurcation points at the single-mode response curve composed of only u 3 . This leads to two equations below, and Similarly, combining Eqs. 35 and 36 to find the position and number of their intersection points, we can predict the occurrence of modal interaction for the only mode 3 forced situation. Figure 7 shows the position relationship between these two curves of the  Fig. 6 The curves of Eqs. 33 (red) and 34 (green) for three different magnitudes of force amplitude cases: a P m2 = 0.4 × 10 −3 , b P m2 = 0.9 × 10 −3 and c P m2 = 1.5 × 10 −3 . The system parameters are: ω n2 = 1.005, ω n3 = 1.015, μ = −0.05, and ζ ≈ 0.001. The red dots mark the intersection points and the black solid and dashed curves represent the stable and unstable force response curves generated using AUTO. (Color figure online) softening system for two force amplitude situations: P m3 = 0.5 × 10 −3 and P m3 = 1.0 × 10 −3 . The intersection points are marked by red dots. The corresponding backbone curves S3 and S5 ± in Fig. 4 and response curves computed by AUTO are also shown to make the comparison.
For the two forced situations, there are zero intersection point, in Fig. 7a, and two intersection points, in Fig. 7b, between curves of Eqs. 35 and 36, respectively. As expected, there is no mixed-mode response in Fig. 7a. In Fig. 7b, it can be seen that the mixedmode response curves following S5 ± emanate from the single-mode response curve from the intersection points. In fact, as there is only one intersection point between the backbone curve S3 and the curve of Eq. 36, the single-mode response curve of u 3 , Eq. 35, can intersect with the curve of Eq. 36 at two points at most in the interesting frequency range. This means that only one kind of mixed-mode resonant responses which is the one following S5 ± shown in the results may bifurcate from the single-mode response curve of u 3 when the linear mode 3 of the system is excited. This is also in agreement with the expectation of the backbone curve results that only the out-of-unison resonant backbone curves S5 ± emanate from S3.

Conclusion
In this paper, we have considered the N −1 modal interactions that occur in a three-degree-freedom lumped mass system. In particular we considered the potential modal interactions of the system by analysing the backbone curves of the undamped, unforced system. This is an important topic because in lightly damped structures the dynamic behaviour is largely determined by the properties of the underlying undamped dynamic system.
First the undamped, unforced case was considered. In particular the modal interaction case that occurs when all the underlying linear modal frequencies are close was examined (i.e. ω n1 : ω n2 : ω n3 ≈ 1 : 1 : 1). In this case the first mode is linear because of the symmetry of the system and the other two modes will potentially interact with each other when the special parameters are chosen. We showed how this system can be analysed using a normal form transformation to obtain the nonlinear backbone curves of the undamped, unforced response. Following this, the response in the  Fig. 7 The curves of Eq. 35 (red) and Eq. 36 (green) for two different magnitudes of force amplitude cases: a P m3 = 0.5 × 10 −3 and b P m3 = 1.0 × 10 −3 . The system parameters are: ω n2 = 1.005, ω n3 = 1.015, μ = −0.05, and ζ ≈ 0.002. The red dots mark the intersection points. The grey curves represent the corresponding backbone curves S3 and S5 ± and black solid and dashed curves represent the stable and unstable forced response curves generated using AUTO. (Color figure online) frequency domain of the corresponding lightly damped and harmonically forced system was obtained using the continuation software AUTO-07p. This results were compared with the backbone curves to show its validity for predicting the nonlinear resonant behaviour of the system during N − 1 modal interactions.