Subscribe Now Subscribe Today
Research Article

Lyapunov Exponents Estimation of the Impact Oscillation Using Observed Time Series

June-Yule Lee
Facebook Twitter Digg Reddit Linkedin StumbleUpon E-mail

This study concerns the estimation of the Lyapunov exponents for impact oscillation using observed time series. An equivalent dynamical system is constructed by delay embedding from impacting series and the averaged exponents is computed by applying QR decomposition. Finally, the results show an agreement with the previous theoretical results.

Related Articles in ASCI
Similar Articles in this Journal
Search in Google Scholar
View Citation
Report Citation

  How to cite this article:

June-Yule Lee , 2006. Lyapunov Exponents Estimation of the Impact Oscillation Using Observed Time Series. Journal of Applied Sciences, 6: 1604-1607.

DOI: 10.3923/jas.2006.1604.1607



Impact dynamics is considered to be one of the most important problems in mechanical vibrating systems. Such impacting oscillation may occur in the motion with amplitude constraining stops (Fig. 1). This problem has been researched in detail using mathematical model (Shaw and Holmes, 1983; Budd et al., 1995; Luo and Xie, 2001) and the simulation results showed in Fig. 2 (Lee, 2005). Different types of impacting response can be predicted from bifurcation diagram for varying the control parameter. For the purpose of condition monitoring of the system, Lee (2006) calculated the spectrum of Lyapunov exponents for varying the driving frequency. The results showed an agreement with the bifurcation diagram of Fig. 2. In this study, the estimation of Lyapunov exponents using the impacting time series is considered, where the delay embedding techniques proposed by Darbyshire and Broomhead (1996) and Babrook et al. (1997) are applied. In the following study, a complete procedure of computing Lyapunov exponents based on time series is developed and the results are compared to the previous theoretical results.


An impact oscillator is considered in Fig. 1, where the system is under a harmonic excitation with an amplitude constraint. The equation of motion, when not impacting, is


where are the acceleration, the velocity, the displacement, the mass, the damping, the stiffness, the forcing amplitude, the forcing frequency and the amplitude constraint, respectively.

Fig. 1: Impact oscillator

Fig. 2:
Bifurcation diagram of the impact oscillator for the control parameters m = 1, c= 0.01, k = 1, p = 1, r = 0.8, g = 0 and ω = 0 to 5

The impact occurs whenever x = g and the velocity is modeled as v(t+) = -rv(t-) where r is restitution coefficient, (t+) is the time after impact and (t-) is the time before impact.

Fig. 3:
The impact oscillator for the control parameters m = 1, c = 0.01, k = 1, p = 1, r = 0.8, g = 0 and ω = 2 (a) Displacement (b) Phase trajectory (c) Impact map (d) Estimated Lyapunov exponents using the observed time series of impact velocity, where the nearby vectors Nb=5 and the embedding dimension m = 2

The result shows in Fig. 2 is a bifurcation diagram with driving frequency range from 0 to 5. The regions of the stable and unstable impact oscillations are clearly demonstrated. For the case of the driving frequency at ω = 2, the displacement of the impact oscillation is shown in Fig. 3a. This is an example of a periodic impact oscillation and the phase trajectory demonstrates a close orbit in Fig. 3b. For the same control parameters but at driving frequency ω=2.8, the result shows in Fig. 4a and b demonstrate a non-periodic impact oscillation.


The Lyapunov exponents are measures of the average rate of the convergence or divergence of the system and defined as


This is based on m-dimension system with the mapping I:Rm → Rm,


For a small perturbation of this mapping tn ⇒ tn+δtn, the Taylor expansion of Eq. 3 and the linearization model is


The stability at point tn is decided by the eigenvalues of the Jacobian matrix DI(n)(tn). Thus m exponents of the average rate λm can be calculated if the Jacobian matrix is available.

In the following investigation, the estimation of the Jacobian matrix is considered using the simulated impacting time series. Based on the definition of the Lyapunov exponents in Eq. 2, the exponents for impact oscillation are


where DI(n)(tn) is the Jacobian matrix for impact map


Fig. 4:
The impact oscillator for the control parameters m = 1, c = 0.01, k = 1, p = 1, r = 0.8, g = 0 and ω = 2.8 (a) Displacement (b) Phase trajectory (c) Impact map (d) Estimated Lyapunov exponents using the observed time series of impact velocity, where the nearby vectors Nb = 5 and the embedding dimension m = 2

Now consider a practical situation that the Jacobian matrix of the impact system is unknown and only a time series of impacting signals is known. Thus, the aim is to estimate Jacobian matrix using observed time series.

For m-dimensional system, an equivalent dynamical system can be constructed using delay embedding from time series (Darbyshire and Broomhead, 1996; Babrook et al., 1997)


and the embedding representation of Eq. 7 using the delay coordinate embedding techniques is


or in matrix form


Thus the estimated Jacobian matrix DIn(tn) of Eq. 6 is accomplished by estimation the difference vector δtn and δtn+1 in Eq. 4. For example, DI(t1) is calculated by examining how a local area evolves from one vector t = [t1t2t3…tm] to another t = [t2t3t4…tm+1] in Eq. 9, which is achieved by finding several nearby vectors Nb and forming a neighborhood matrix B1 (Nbxm) and its evolved matrix B2 (Nbxm) . Thus the Eq. 4 is expressed by


and the pseudo-inverse solution of the estimated Jacobian matrix is given by


where ( )T is transposition matrix and ( )-1 is inverse matrix.

In computation of the eigenvalues of the matrix J in Eq. 5, that is given by


Only a few iterations, the matrix J becomes very large for chaotic case and null for the periodic case because the product of the matrix DI in Eq. 12. To overcome this problem, the QR factorisation techniques in proving the multiplicative ergodic theory (Eckmann and Rulle, 1985; Johnson et al., 1987) are applied. Given a matrix A, there is a unique factorization


where R is a square right upper triangular matrix and Q has orthogonal columns. Consider an arbitrary orthogonal matrix Q0 and apply the orthogonalisation procedure to the impact map gives


The columns Q1 are now an orthogonal basis for the impact map at n = 1 and R0 gives the relationship between this map and the image of the Q0 under the Jacobian matrix DI(1). This process can be repeated recursively and gives


This method removes the need to compute products of large numbers of matrices. Thus the Lyapunov exponents have the form of a simple ergodic average and can be expressed by


Summary of the processes of estimating Lyapunov exponents from time series:

1. Embed the data using Eq. 9.
2. Locate nearest embedded vector and form a neighbourhood matrix B2.
3. Evolve one step and locate nearest embedded vector and form a neighbourhood matrix B2.
4. Calculate the estimated Jacobian matrix using the pseudo-inverse method in Eq. 11.
5. Perform QR decomposition.
6. Calculate the average exponents using Eq. 16.
7. Repeat the processes (2-6) until the exponents converge to fixed values.


The result shows in Fig. 3c is the points set in impact map at driving frequency ω = 2 and the estimated Lyapunov exponents are shown in Fig. 3d. The Lyapunov exponents converge to a non-positive value after 5000 iterations. This indicates that the impact oscillation is stable and is in consistent with Fig. 2. For the case of driving frequency at ω=2.8, the results show in Fig. 4c. The “strange” point sets in impact map demonstrate the impact oscillator is chaos. One positive Lyapunov exponents shown in Fig. 4d confirms that the impact oscillation is unstable. Both results in Fig. 3d and 4d are also in consistent with Fig. 2.


The algorithms of computing Lyapunov exponents were investigated using observed impacting series. The calculating processes were (1) estimated Jacobian matrix by embedding techniques and (2) averaged the exponents by applying QR decomposition. Finally, the results of the estimated Lyapunov exponents showed an agreement with the bifurcation diagram.

1:  Banbrook, M., G. Ushaw and S. McLaughlin, 1997. How to extract Lyapunov exponents from short and noisy time series. IEEE Signal Process., 45: 1378-1382.
CrossRef  |  

2:  Budd, C., F. Dux and A. Cliffe, 1995. The effect of frequency and clearance variations on single-degree of freedom impact oscillators. J. Sound Vib., 184: 475-502.
Direct Link  |  

3:  Darbyshire, A.G. and D.S. Broomhead, 1996. Robust estimation of tangent maps and Lyapunov spectra. Phys. D. Nonlinear Phenomena, 89: 287-305.
CrossRef  |  

4:  Eckmann, J.P. and D. Rulle, 1985. Ergodic theory of chaos and strange attractors. Rev. Modern Phys., 57: 617-656.
CrossRef  |  

5:  Johnson, R.A., K.J. Palmer and G.R. Sell, 1987. Ergodic properties of linear dynamical systems. SIAM J. Math. Anal., 18: 1-33.
CrossRef  |  

6:  Lee, J.Y., 2005. Motion behavior of impact oscillator. J. Mar. Sci. Technol., 13: 89-96.
Direct Link  |  

7:  Lee, J.Y., 2006. The lyapunov exponents of the impact oscillator. J. Applied Sci., 6: 80-84.
CrossRef  |  Direct Link  |  

8:  Luo, G. and J. Xie, 2001. Bifurcation and chaos in a system with impacts. Physica D: Nonlinear Phenomena, 148: 183-200.
CrossRef  |  

9:  Shaw, S.W. and P.J. Holmes, 1983. A periodically forced linear piecewise linear oscillator. J. Sound Vib., 90: 129-155.
CrossRef  |  

©  2021 Science Alert. All Rights Reserved