Modelling of nonlinear modal interactions in the transient dynamics of an elastic rod with an essentially nonlinear attachment

Abstract We perform system identification and modelling of the strongly nonlinear modal interactions in a system composed of a linear elastic rod with an essentially nonlinear attachment at its end. Our method is based on slow/fast decomposition of the transient dynamics of the system, combined with empirical mode decomposition (EMD) and Hilbert transforms. The derived reduced order models (ROMs) are in the form of sets of uncoupled linear oscillators (termed intrinsic modal oscillators – IMOs), each corresponding to a basic frequency of the dynamical interaction and forced by transient excitations that represent the nonlinear modal interactions between the rod and the attachment at each of these basic frequencies. A main advantage of our proposed technique is that it is nonparametric and multi-scale, so it is applicable to a broad range of linear as well as nonlinear dynamical systems. Moreover, it is computationally tractable and conceptually meaningful, and it leads to reduced order models of rather simple form that fully capture the basic strongly nonlinear resonant interactions between the subsystems of the problem.


Introduction
The need for system identification and reduced order modeling arises from the fact that, presented with sensor data, the analyst is generally unaware of details of the underlying dynamical system from which they originated. The common approach to this dilemma is to assume that the dynamical system is linear and that the measured responses are stationary in time. This facilitates the use of the numerical Fourier transform (FT) followed by experimental modal analysis (EMA) [4] to extract natural frequencies, mode shapes and modal damping ratios, from which the parameters of the assumed linear model can be determined once the mass distribution is known. This approach, which is fully nonparametric, has served the dynamics and controls community well, even in the presence of weakly nonlinear system behavior. Clearly, though, as systems become more complex, lightweight and flexible, incorporating not only electrical and mechanical components but also additional elements and requirements such as, for example, bio-inspired features, the likelihood exists that the underlying behavior will be strongly nonlinear and nonstationary. Such dynamical characteristics can result from, for example, local 1 buckling, plastic deformations, clearance and backlash, hysteresis, friction, and so forth, and as we think more in terms of multi-physics problems, one must also include nonlinearities due to interfacial effects such as shear lag between actuator and structure, fluid-structure interactions, sensor-tissue interactions, and so forth. In general, a physics -based parametric model of the system will not be known a priori.
However, given a sufficiently dense set of sensors, measured time series recorded throughout the system will contain all of the information reflecting both nonlinearity and nonstationarity. The classical FT is not able to properly isolate and extract this information and, in fact, may lead the less experienced analyst to mistake phenomena such as internal and combination resonances for natural frequencies, to fail to account for sensitivities of the response to force and voltage magnitudes, initial conditions, and to miss or misinterpret other unique behaviors. These observations highlight the importance of developing effective, straightforward, nonparametric system identification and reduced order modeling methods for characterizing strongly nonlinear, complex, multi-component systems that will be as utilitarian as (the well established) EMA is for linear systems. The need for developing such a nonlinear system identification technique of broad applicability is dictated by the limitations of current system identification techniques which are either applicable to only linear systems, or are tailored to special classes of smooth nonlinear systems [11]. The difficulty in developing nonlinear system identification methodologies valid for broad classes of dynamical systems is due to the well-recognized highly individualistic nature of nonlinear systems, a feature that restricts the unifying dynamical features which are amenable to system identification [16].
Our proposed nonparametric nonlinear system identification (NNSI) technique relies solely on direct time series measurement and post-processing, and leads to nonlinear interaction models (NIMs) of simple form, represented as sets of uncoupled linear oscillators with time invariant coefficients forced by inhomogeneous terms representing nonlinear modal interactions. Key to our method is a slow/fast partition of time series measurements that leads to the identification of the basic (dominant) fast frequencies of the dynamics (which also govern the dimensionality of the resulting reduced order model), and slow flows governing the important (essential) dynamics of the problem. We then perform direct system identification of the slow flow dynamics employing advanced post-processing computational algorithms, namely, the wavelet x k ðtÞ instantaneous frequency transform, the Hilbert transform, and empirical mode decomposition (EMD). In particular, EMD is a signal analysis method which has received much attention lately due to its application in a number of fields [7,8], but up to now was regarded as lacking a theoretical or formal foundation. As shown, however, in [14,15] a solid theoretical foundation for EMD can be gained in terms of the slow flow dynamics of the problem. Based on this key finding our NNSI technique successfully integrates the aforementioned analytical and computational algorithms in a powerful and robust synergistic algorithm for constructing ROMs in broad classes of nonlinear dynamical systems.

Transient resonance captures (TRCs) and nonlinear damped transitions
We consider an unforced viscously damped dispersive finite rod on a distributed elastic support, possessing an essentially nonlinear attachment at its right end (cf. Fig. 1). The essential nonlinearity of the attachment is due to the lack of a linear component in the stiffness connecting it to the rod and, as discussed below, leads to strongly nonlinear modal interactions between the rod and the attachment. The attachment consists of a light mass (scaled by the small parameter 0 < e ( 1Þ and is connected to the rod by the essentially nonlinear stiffness which is in parallel to a viscous damper. Denoting by vðtÞ and uðx; tÞ the responses of the nonlinear attachment and the rod, respectively, we express the equations of motion of the system in the form @ 2 uðx; tÞ @t 2 þ x 2 0 uðx; tÞ þ ek 1 @uðx; tÞ @t À @ 2 uðx; tÞ @x 2 ¼ 0; 0 6 x 6 L uð0; tÞ ¼ 0; @uðL; tÞ @x ¼ Àe€ vðtÞ e€ vðtÞ þ ek 2 _ vðtÞ À @uðL; tÞ @t ! þ C½vðtÞ À uðL; tÞ 3 ¼ 0 where rðxÞ and sðxÞ are the initial displacement and velocity distributions of the rod, and the viscous dissipative terms of the system are assumed to be small of OðeÞ. The normalized frequency x 0 represents the cut-off frequency of the uncoupled dispersive rod of infinite extent, i.e., the bounding frequency separating travelling waves from near-field standing wave solutions in the rod. The strongly nonlinear transient dynamics of this system was computationally studied in previous work [5,20,6,18,21] by constructing the frequency -energy plot (FEP) of periodic orbits of the undamped Hamiltonian system. Damped transitions in the dynamics of the system can then be studied by superimposing their wavelet transform spectra on the Hamiltonian FEP, recognizing that for sufficiently weak damping the effect of damping is purely parasitic, its only function being the determination of the branches of periodic or quasiperiodic solutions in the FEP that are 'visited' by the dynamics. Due to the essential stiffness nonlinearity of the attachment, system (1) may exhibit very complicated transient dynamics indeed, including multi-frequency transitions between different response regimes.
In this work we will focus on one such specific damped transition. For this, the transient damped dynamics of system (1) was computed by the finite element method (employing a 200 element discretization), under the assumption of proportional viscous damping distribution (with a viscous damping matrix of the form D ¼ a 1 M -where M is the mass matrix of the discretized model of the rod), and parameters e ¼ 0:05; C ¼ 1:0; L ¼ 1:0; x 0 ¼ 1:0; k 2 ¼ 0:2; a 1 ¼ 0:005. The leading eigenfrequencies and viscous damping ratios of the linear dispersive rod with no attachment are listed in Table 1. The specific initial conditions for the damped response analyzed in this work are given by 1=5 (cf. Fig. 2c). Following [20,18], a tongue labeled as T ðnÞ p=q denotes the branch of subharmonic motions of the system where the dominant harmonic component of the nonlinear attachment is nearly equal to ðp=qÞx n , whereas that of the rod is x n , where x n is its n-th linearized eigenfrequency; it follows that on the subharmonic tongue T ðnÞ p=q the periodic response of the system of Fig. 1 possesses two main harmonics at frequencies x n and ðp=qÞx n .
Following this notation, and referring to Fig. 2c, the damped transient motion considered in this study initially possesses two dominant harmonic components at frequencies x 4 (the dominant harmonic of the rod) and x 4 =5 (the dominant harmonic of the nonlinear attachment), but as energy decreases due to damping dissipation, multi-frequency transitions occur, as evidenced by the irregular amplitude modulations of the corresponding time series (especially the one of the nonlinear attachment -cf. Fig. 2b). These multi-frequency transitions are studied in Fig. 2c, where the wavelet spectrum of the relative motion ½vðtÞ À uðL; tÞ is superimposed on the Hamiltonian FEP of system (1). In this case there exist three distinct stages of the damped motion. In the initial stage, there is a high -energy transition from the subharmonic tongue T ð4Þ 1=5 (where the motion of the system is initiated) to tongue T ð1Þ 1=3 . This is followed by prolonged subharmonic transient resonance capture (TRC) of the dynamics on T ð1Þ 1=3 [21] as signified by the presence of strong harmonics in the response at frequencies x 1 =3 and x 1 . Finally, there is a low-energy transition to a linearized state, where the response of the nonlinear attachment is negligible, and the dynamics is dominated by the response of the linear rod.
We wish to perform system identification and reduced order modeling of the described strongly nonlinear modal interactions between the rod and the nonlinear attachment. To address this task we will rely solely on direct analysis of the time series of the rod end and the nonlinear attachment. In addition to the wavelet transform we will employ additional advanced post-processing tools; namely, empirical mode decomposition (EMD) and Hilbert transform analysis. In the following sections we will discuss the basic steps of our analysis in detail.

Empirical mode decomposition (EMD), Hilbert spectrum analysis and masking signals
Since the empirical mode decomposition (EMD) method combined with Hilbert spectral analysis was introduced in [7,8], numerous applications to system identification [23,17], nonlinear dynamics [20,5,18,6] and damage detection [24,10] have appeared. This decomposition method, based on identifying the characteristic time scales in measured time series, is adaptive, highly efficient, and suitable for nonlinear and nonstationary processes. In particular, EMD yields a complete and nearly (but not completely) orthogonal basis of intrinsic mode functions -IMFs; these are oscillatory modes embedded in the time series, each with its characteristic time scale.
Until now EMD was applied for time series analysis in an ad hoc fashion. The main loop of the EMD algorithm for extracting the IMFs from a signal xðtÞ consists of the following steps [7,8]: (i) identify all extrema of xðtÞ; (ii) perform (spline) interpolations of the minima and maxima of xðtÞ, ending up with two envelopes e minðtÞ and e maxðtÞ, respectively; (iii) compute the average curve RðtÞ ¼ e minðtÞ þ e maxðtÞ ½ =2 (as a residual); (iv) extract the remainder cðtÞ ¼ xðtÞ À RðtÞ; (v) apply the previous algorithm on the reminder cðtÞ until the residual RðtÞ can be considered as zero-mean under some tolerance (i.e., a stopping criterion). Once this criterion (through the sifting process) is met, the remainder cðtÞ is regarded as the effective IMF. By subtracting this IMF from the original time series and applying the algorithm iteratively we extract additional IMFs, so that the original signal xðtÞ is decomposed sequentially from high-to low-frequency components as: By this construction process the superposition of the K leading IMFs reconstruct approximately the measured time series; however, due to the ad hoc nature of the sifting algorithm only a subset of these IMFs are physically meaningful with the rest being of spurious nature. As discussed in [20,5,6], the dominant (and physically meaningful) IMFs can be identified by comparing their instantaneous frequencies to the wavelet spectra of the original time series: the instantaneous frequencies of the dominant IMFs coincide with the dominant harmonics of the wavelet spectra. This process also identifies the dominant time scales (or frequencies) of the dynamics in the time series. The instantaneous frequency of an IMF is computed by complexification and application of the Hilbert transform. The Hilbert transform (HT) hðtÞ of a (mono-component) signal yðtÞ is computed as where PV stands for Cauchy principal value. Moreover, a complex function whose imaginary part is the Hilbert transform of the real part is analytic. Motivated by this result we may complexify the k-th IMF c k ðtÞ of the time series xðtÞ by defining the analytic complex function w k ðtÞ c k ðtÞ þ jH½c k ðtÞ ð5Þ where j ¼ ðÀ1Þ 1=2 . This enables the computation of the instantaneous amplitude and phase of the k-th IMF as from which the instantaneous frequency of the IMF is computed asx k ðtÞ ¼ _ h k ðtÞ. Although these relations define the instantaneous frequency of the IMF regardless of the bandwidth, it has been observed that the notion of instantaneous frequency has physical meaning only for narrowband signals; e.g., for high-frequency modulated signals. However, it is possible to implement frequency demodulation for wideband signals as well, so that the instantaneous frequency can be considered as an average of all frequencies that exist in the IMF at a given instant, while the instantaneous bandwidth can be considered as the deviation from that average [19]. As discussed in [14], application of EMD to time series with dominant frequencies whose ratios are in the range from 0.5 to 2 Hz lead to 'complicated-mode' IMFs exhibiting beat phenomena; the issue of how to decompose such signals with closely spaced frequencies is important and not yet completely resolved.
Although EMD combined with the Hilbert transform (HT) forms a powerful post-processing tool for extracting intrinsic oscillating components from a time series and identifying the dominant time scales of the dynamics, it has some important deficiencies. First, application of EMD may lead to spurious IMFs, so that physically meaningful results from EMD can only be obtained if these are omitted from further consideration in the analysis; spurious IMFs are the direct result of the well established lack of orthogonality of the IMFs. The deletion of spurious IMFs from further consideration can be performed by employing WT spectra as discussed previously. Second, there are concerns regarding the frequency resolution of the EMD results. Indeed, in order to obtain meaningful results when applying HT to the IMFs, it is necessary that these are mono-component or, at least, narrowband (otherwise mixed-mode IMFs in the form of beat phenomena are obtained possessing closely spaced frequencies). Finally, there are issues concerning the uniqueness of the EMD results. EMD does not result in a unique decomposition of a measured time series since it is applied in an ad hoc manner and depends on a free stopping parameter; that is, EMD is not robust in practice. The set of extracted IMFs can be considered as a basis for reconstructing the original measured time series if it satisfies (or nearly satisfies) the basic conditions of completeness and orthogonality. By virtue of the EMD sifting algorithm, completeness of the IMFs is guaranteed by construction. It is the lack of orthogonality between IMFs, however, that generates spurious features in the results and prevents uniqueness of the decomposition.
The lack of IMF orthogonality, however, is a motivation for devising ways in order to refine EMD. First, we note that orthogonality of the IMFs resulting from EMD analysis can be quantified; e.g., in (3) we may neglect the residual term R Kþ1 ðtÞ and introduce the square of the truncated expression (3) Then, an index of orthogonality (IO) can be defined by employing the cross term in (7) so in the limit of complete orthogonality of the IMFs it holds that IO ! 0. This observation will help us refine the EMD process by requiring that the derived IMFs correspond to the least value of IO [typically below a prescribed tolerance, e.g., of Oð10 À3 Þ]. This will lead to the development of an advanced EMD process -AEMD [14] through the use of masking signals [for alternative techniques for improving EMD refer to [1,2,9,12,22]].
Deering and Kaiser [3] claimed and demonstrated that utilizing masking signals can improve EMD's deficiency in frequency resolution, which should be useful for intermittency (or mode mixing). The main idea of using a masking signal on the measured time series is to offset the inability of standard EMD to generate proper IMFs (that is, nearly orthogonal and narrowband monocomponent), so the success of the AEMD depends on the proper choice of the amplitude and frequency of the applied masking signals. In addition, application of masking signals can successfully address the existence of inflection points in the time series which cannot be accounted for in the standard EMD [14].
Application of the AEMD technique will be demonstrated with an example which is reproduced from [14] that illustrates its main features. We consider the signal xðtÞ ¼ cosð2pf 1 tÞ þ 0:3 cosð2pf 2 tÞ, f 1 ¼ 1 Hz, f 2 ¼ 3 Hz, for which the standard EMD cannot separate the weak higher harmonic component from the strong lower frequency one. A schematic explaining the application of the masking signal xmaskðtÞ ¼ 0:4x max cosð2pf max tÞ (where f max ¼ 20f 1 and x max is the maximum amplitude of xðtÞÞ is depicted in Figs. 3 and 4. To perform AEMD of xðtÞ, we perform standard EMD analysis of the signal x þ ðtÞ xðtÞ þ xmaskðtÞ yielding the IMFs c þ k ðtÞ, k ¼ 0; 1; 2 (cf. Fig. 3a). Then we perform standard EMD analysis of the signal x À ðtÞ xðtÞ À xmaskðtÞ yielding the corresponding set of IMFs c À k ðtÞ, k ¼ 0; 1; 2 (cf. Fig. 3b). Finally, we compute the mean curves of the respective IMFs to minimize the influence of the masking signal; i.e., we define the average IMFs Apparently, the leading IMFs c AE 0 ðtÞ correspond to the harmonic of the masking signal at frequency f 3 ¼ 20f 1 , which are out-of-phase with each other; therefore, averaging this set of IMFs leads to c 0 % 0, which implies that the effects of the masking signal are almost completely eliminated from further analysis (cf. Fig. 4). The other two IMFs c 1 and c 2 provide a physically meaningful decomposition in terms of nearly monochromatic components at frequencies f 1 and f 2 , with small index of orthogonality equal to IO ¼ 4:7 Â 10 À3 .
We note that, the selection of the masking signal (or of its parameters) is not unique; e.g., using different values for the frequency f 3 in the previous example may also lead to physically meaningful AEMD results. Examination of several masking signals to perform AEMD analysis provides some general guidelines for utilizing masking signals in practical signal analysis. First, the frequency of the masking signal need not necessarily be high. The presence of higher frequency components in the response implies smaller time step sizes and, hence, larger number of data points required for signal processing, which causes long computational time and sometimes ill-behaved decomposition. In fact, the use of a masking signal has the effect Fig. 3. Application of a masking signal xmaskðtÞ ¼ 0:4xmax cosð2pfmaxtÞ; f max ¼ 20f 1 to the signal xðtÞ ¼ cosð2pf1tÞ þ 0:3 cosð2pf2tÞ; f 1 ¼ 1 Hz; f 2 ¼ 3 Hz: EMD analysis of (a) xðtÞ þ xmaskðtÞ, (b) xðtÞ À xmaskðtÞ [14]. of adding artificial noise to the original time series, magnifying the hidden (or weak) frequency components. Therefore, the choice of the masking signal should hinge on whether it contains frequency components which are similar to those of the original signal to be analyzed. This requires a pre-test (such as Fourier or wavelet spectral analysis) of the original signal to identify its frequency content. Second, a successful AEMD of a signal clearly depends on the masking signal amplitude; that is, it is wise to select the amplitude of the masking amplitude to be sufficiently high in order to amplify sufficiently hidden or weak frequency components in the original time series.
We analyze now the transient multi-frequency damped transition of system (1) depicted in Fig. 2. In Fig. 5 we present the power spectrum and corresponding wavelet transform of the response of the rod end (cf. Fig. 2a), and notice that it possesses two main distinct harmonics. Hence, it is not necessary to apply a masking signal to perform EMD in this case, so the results of application of standard EMD analysis are depicted in Fig. 6. There are two dominant, mono-component IMFs with index of orthogonality IO ¼ 0:0031. These are almost proper IMFs (according to the discussion above), and their superposition accu-  8 rately reconstructs the rod end response (cf. Fig. 6c). Moreover, the instantaneous frequencies of the dominant IMFs coincide with the dominant harmonics of the wavelet transform spectra of the rod end response (cf. Fig. 6a and 6b). The two dominant IMFs possess nearly constant frequencies identical to the first and fourth eigenfrequencies of the rod with no attachment. This correlates fully with the representation of the wavelet spectrum of this damped transition on the FEP (cf. Fig. 2c). Studying the power spectrum and corresponding wavelet transform of the response of the nonlinear attachment depicted in Fig. 7, we note that there exist closely spaced, time-varying harmonics in this case, so we need to apply a masking signal in order to obtain nearly mono-component and nearly orthogonal (i.e., proper) IMFs. In this case we perform AEMD by applying the masking signal xmaskðtÞ ¼ À0:4V max cosð2pf 1mask tÞ À 0:36V max cosð2pf 2mask tÞ to the time series, where V max is the maximum response of the response of the nonlinear attachment, 2pf 1mask ¼ 2:214 rad=s % x 4 =5, and f 2mask ¼ 4f 1mask . In Fig. 8 we depict the resulting (almost proper) IMFs together with their instantaneous frequencies; the index of orthogonality of these IMFs is equal to IO ¼ 0:0139. Similar to the response of the end of the rod, the response of the nonlinear attachment possesses two dominant IMFs that accurately reconstruct the original time series (cf. Fig. 8c and 8f). Considering, however, the instantaneous frequencies of the two IMFs, we note that in this case the first dominant IMF possesses a nearly constant frequency equal to one-fifth of the fourth eigenfrequency of the rod with no attachment, whereas the second IMF possesses a slowly varying eigenfrequency which ends up 'locking' to one third of the first eigenfrequency of the rod in the time interval 250 < t < 500.
These results correlate fully with the damped transitions depicted in the FEP of Fig. 2c, as they indicate that the first IMF of the rod end response and the first IMF of the nonlinear attachment capture the initial 1:5 TRC of the damped dynamics during initiation of the motion on the subharmonic tongue T 1=3 (with the rod oscillating at frequency x 1 and the nonlinear attachment at frequency x 1 =3). Clearly, there are additional (secondary) TRCs and transitions in the dynamics, especially in the initial high-energy regime of the motion, but in our analysis we will focus only on the aforementioned two dominant TRCs.
Hence, an interesting decomposition of the two dominant TRCs of the dynamics in terms of the dominant IMFs of the rod and attachment responses is achieved by AEMD. In the next section we demonstrate how the AEMD results can be employed to construct an accurate reduced order model of the strongly nonlinear damped transition considered above.

Nonlinear nonparametric system identification (NNSI)
Based on the previous AEMD results we perform nonparametric nonlinear system identification leading to nonlinear interaction model (NIM) for the rod-attachment nonlinear damped transition of Fig. 2. Considering first the rod end response and referring to the AEMD results depicted in Fig. 6, we note that this response can be modeled by a system of uncoupled, forced linear oscillators, termed intrinsic modal oscillators -IMOs where cc denotes complex conjugate. The inhomogenous (forcing) terms in (9) are in the form of 'fast' oscillating terms e jx k t , k ¼ 1; 4 modulated by the 'slowly' varying complex amplitudes K u k ðtÞ, k ¼ 1; 4; this implies that we have assumed that the modulations K u k ðtÞ vary much slower than the corresponding carrying signals e jx k t . The main motivation for considering this form of ROM for the rod end response, is so that the response of each of the two IMOs reproduces approximately one of the dominant IMFs of the dynamics of the rod end (cf. Fig. 6). Since the superposition of the IMFs reconstructs the original time series, the same result should hold for the combined response u 1 ðtÞ þ u 2 ðtÞ of the IMOs in (9). The reasoning behind the specific structure of the ROM (9), i.e., in the form of a set of linear damped and forced oscillators, lies in the fact that each of the dominant IMFs of the rod end response possesses two constant dominant fast frequencies equal to the first and fourth eigenfrequencies of the rod, respectively. It follows that (at least in principle) each of the dominant IMFs can be regarded as the response of a damped linear oscillator (i.e., an IMO) with an eigenfrequency equal to one of the fast frequencies of the rod end response and a viscous damping coefficient k u k , k ¼ 1; 2. The specific forms of the inhomogeneous (forcing) terms in (9), i.e., as modulated periodic signals, are dictated by the realization that any alternative type of excitation would be off-resonance and, hence, their effects on the rod dynamics would be small, of a secondary nature. In essence, these inhomogeneous terms represent the nonlinear modal interactions between the rod end and the nonlinear attachment at the dominant fast frequencies of the damped response of the rod edge. We also observe that the IMOs of the ROM (9) are uncoupled since the dominant fast frequencies of the rod end response are distinct (cf. the wavelet spectra of Fig. 5); when the dynamics possesses closely spaced frequencies, this should be accounted for in the ROM, and the simple form (9) should be adjusted to reflect 'mode mixing' (beat phenomena). Finally, we note that this system identification and reduced order modeling process is fully nonparametric in the sense that it does not depend on the specific form of nonlin-  . AEMD analysis of the response of the nonlinear attachment, 0 < t < 1000: (a) 1st IMF and its instantaneous frequency superimposed on the WT spectrum, (b) 2nd IMF and its instantaneous frequency superimposed on the WT spectrum, (c) reconstruction of response by superimposing the two leading (dominant) IMFs, (d) detail of 1st IMF and its instantaneous frequency for 0 < t < 100, (e) detail of 2nd IMF for 0 < t < 400 and its instantaneous frequency for 0 < t < 300, (f) detail of reconstruction of the response by the two dominant IMFs for 0 < t < 200. earity of the system, as it relies on direct analysis of measured time series. Moreover, the resulting NIM is of rather simple form and the strongly nonlinear modal interactions are represented in the most general form, i.e., as time dependent functions. Following similar reasoning and taking into account that the damped response of the attachment possesses two dominant fast frequencies approximately equal to x 4 =5 and x 1 =3, we construct the approximate ROM for the damped dynamics of the nonlinear attachment Each of the two IMOs of the ROM (10) reproduces approximately a dominant IMF of the response of the nonlinear attachment. From the AEMD results of the previous section we recall that the first IMF possesses a nearly constant frequency equal to one-fifth of the fourth eigenfrequency of the rod, whereas the second IMF possesses a slowly varying eigenfrequency that 'locks' to one third of the first eigenfrequency of the rod at a later-stage of the motion. In the ROM (10) the time variation of the instantaneous frequency of the second IMF is not taken into account; but rather, to a first order of approximation, the NIM possesses constant eigenfrequencies. This approximation will be reflected in the accuracy of the reproduction of the IMFs by the NIM, especially in the early-and later-stages of the motion of the nonlinear attachment. We note that the first IMO in (10) models the initial stage of the response of the nonlinear attachment, i.e., the initial 1:5 TRC of the damped dynamics during initiation of the motion on the subharmonic tongue T ð4Þ 1=5 ; whereas the second IMO models the second (delayed) TRC of the dynamics on the subharmonic tongue T ð1Þ 1=3 . We make the remark at this point that the basis for the aforementioned reduced order modeling procedure is the construction of NIMs that accurately reproduce the dominant IMFs of the rod end and nonlinear attachment responses. This is not an ad hoc requirement since, as discussed in [14], the slowly varying envelopes and phases of the dominant IMFs of the responses of coupled oscillators coincide with the slow flow of the problem. By slow flow we denote the underlying (governing) dynamics of the problem once the (secondary) fast dynamics are averaged out. Although we did not derive the slow flow of the dynamics of the system considered herein, the responses of the NIMs (9) and (10) approximately reproduce the dominant IMFs, so the slowly varying envelopes and phases of these responses represent the slow flow dynamics governing the rod-attachment nonlinear interaction. Hence, the construction of the NIMs (9) and (10) is physics-based and formal; in fact, the aforementioned correspondence between EMD and slow flow dynamics provides a theoretical foundation and gives physical meaning to the dominant IMFs.
Returning now to the NIMs (9) and (10), we compute the damping coefficients and inhomogeneous terms (representing nonlinear modal interactions) by imposing the requirement that each of the IMOs of the ROMs reproduces a dominant IMF of the rod end or the nonlinear attachment responses. To demonstrate this analytical computation we consider an IMO with the general form where KðtÞ is a slow (complex) modulation of the fast periodic signal e jxt . We analyze the dynamics of (11) by performing a slow/fast partition of the dynamics. To this end, we introduce the new complex variable gðtÞ ¼ _ yðtÞ þ jxyðtÞ, and express the real dependent variable and its derivatives as where overbar denotes complex conjugate. Employing (12) the IMO (11) is expressed in the equivalent complex form _ gðtÞ À jxgðtÞ þ ðk=2Þ½gðtÞ þ gðtÞ ¼ KðtÞ: We now assume that the response of (12) can be decomposed in terms of slow and fast components through the relation where the factor jx in (16) accounts for the slightly different definitions of gðtÞ andŵðtÞ. Relations (14) and (16) determine the complex modulation KðtÞ once the damping coefficient k is estimated through an optimization procedure (i.e., by requiring that the response of the IMO best approximates the corresponding IMF). The aforementioned procedure was applied to determine the nonlinear modal interactions K u;v 1;4 ðtÞ forcing the IMOs (9) and (10); for this computation the IMFs depicted in Figs. 6 and 8 were filtered to eliminate high-frequency effects caused by numerical differentiation Lee et al. [14,15]. The results are depicted in Figs. 9-12 for corresponding damping coefficients k u 4 ¼ 0:18; k u 1 ¼ 0:4; k v 4 ¼ 1:0, and k v 1 ¼ 1:2. In essence, the nonlinear modal interactions provide us with a picture of the tran- Fig. 11. Nonlinear modal excitation of the IMO of the response of the attachment at frequency x4=5: (a) magnitude, (b) real part, (c) imaginary part.   (9) and (10) we can reproduce the original strongly nonlinear transient responses, as depicted in Fig. 15. The early-and later-time discrepancies between the responses of the NIM (10) and the nonlinear attachment (cf. Fig. 14c and 14d) are attributed to the time vari-  ations of the dominant frequencies of the nonlinear attachment in these regimes of the dynamics (cf. the wavelet spectrum of Fig. 7b), which was not taken into account the the NIM (10) since fixed eigenfrequencies were assumed for the IMOs In summary, the IMOs (9) and (10) are in the form of sets of linear, uncoupled oscillators, whereas their forcing terms represent the strongly nonlinear modal interactions between the rod and the essentially nonlinear attachment. These forcing terms provide us with information on the level and (slow) temporal dependence of these interactions. The use of linear NIMs for studying strongly nonlinear dynamical interactions was enabled by the accurate identification of the dominant frequencies at which these interactions occurred, and by the equivalence between the IMFs computed by AEMD and the underlying slow flow dynamics which was demonstrated in [14]. This last observation provides a theoretical foundation for the linear NIMs (9) and (10).
Employing this formulation we are able to identify the nonlinear temporal energy exchanges between the rod and the attachment and the precise time scales (frequencies) at which these interactions occur. Hence, the aforementioned procedure leads to multi-scale system identification. In addition, the derived NIMs are of relatively simple form, so they can readily be used for further study of strongly nonlinear dynamical processes and control of these processes. Moreover, the nonparametric way employed for performing the system identification of the dynamics; the fact that it is based on direct analysis of measured time series and the simple linear form of the resulting NIMs enable application of the proposed approach to a broad class of dynamical systems.

Concluding remarks
A main advantage of our proposed technique is that it is nonparametric, eliminating the necessity for a priori assumption of functional forms for stiffness and damping nonlinearities, which might restrict system identification. Hence, at least in principle, it is applicable to a broad range of linear as well as nonlinear dynamical systems, including systems with smooth or non-smooth nonlinearities (such as clearances, vibroimpacts, and dry friction), and strong (even nonlinearizable) or weak nonlinear effects. This is due to the fact that the proposed method directly analyzes the actual measured time series which contain full information of the dynamics and do not rely on computed characteristics of the signals (such as FT analysis). In addition, it is multi-scale and provides directly a measure of the dimensionality of the underlying dynamics (which can be of much smaller order -indeed, many orders of magnitude smaller) than the dimensionality of the underlying computational model. Finally, the method is computationally tractable, conceptually meaningful, and can be used for the construction of accurate low-order models of the dynamics that fully capture the basic resonant interactions between components that give rise to complex and rich dynamical phenomena, such as sudden nonlinear transitions, formation of instabilities, and multifrequency behavior. The presented nonlinear system identification and reduced order modelling technique although nonparametric, was performed under two basic assumptions. First, that the analyzed signals (time series) possessed distinct 'fast' frequencies; this assumption allowed us to derive NIMs in the form of uncoupled linear oscillators (the IMOs) forced by transient nonlinear modal interactions. When the measured time series possesses closely spaced frequencies, the ROMs must be modified to account for beat phenomena between the modes corresponding to these frequencies; even in this case, however, a slow/fast decomposition of the dynamics can be performed so the basic principles of the outlined nonlinear system identification method still apply. The second basic assumption in deriving the NIMs (9) and (10) was that the IMOs possessed constant eigenfrequencies; yet, as shown in the specific application considered in this work, this might not be the case in practice, so that this assumption might introduce errors in the modeling of the dynamics (as indeed was the case in reproducing the early-and later-stage regimes of the response of the nonlinear attachment in Fig. 15). In order to improve the accuracy of the NIMs in such regimes of the motion it is necessary to consider IMOs with modulated eigenfrequencies. Provided that the modulations of the eigenfrequencies are sufficiently slow (compared to the eigenfrequencies themselves), a slow/fast decomposition of the dynamics can still be performed [13]; so the system identification method can be extended in this case as well. This is the focus of current work by the authors.