Journal of Applied Mathematics

Volume 2013 (2013), Article ID 841690, 13 pages

http://dx.doi.org/10.1155/2013/841690

Research Article

Application of the Characteristic Time Expansion Method for Estimating Nonlinear Restoring Forces

Yung-Wei Chen

Department of Marine Engineering, National Taiwan Ocean University, No. 2, Pei-Ning Road, Keelung 20224, Taiwan

Received 17 May 2013; Revised 29 September 2013; Accepted 30 September 2013

Academic Editor: Magdy A. Ezzat

Copyright © 2013 Yung-Wei Chen. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

This paper proposes a characteristic time expansion method (CTEM) for estimating nonlinear restoring forces. Because noisy data and numerical instability are the main causes of numerical developing problems in an inverse field, a polynomial to identify restoring forces is usually adopted to eliminate these problems. However, results of the way doing are undesirable for a high order of polynomial. To overcome this difficulty, the characteristic length (CL) is introduced into the power series, and a natural regularization technique is applied to ensure numerical stability and determine the existence of a solution. As compared to previous solutions presented in other researches, the proposed method is a desirable and accurate solver for the problem of restoring the force in the inverse vibration problems.

1. Introduction

Identification of nonlinear dynamical system is a kind of inverse problems and is usually encountered in engineering applications. For instance, to specify the parameters of dynamical systems is necessary in optimal processes. It is important to analyze and determine the parameters of the system using experimental testing and numerical methods. However, uses of these methods might lead to some challenging problems in the structural mechanic field because a small measurement error can cause large errors in the results of parameters’ identification.

To overcome these inverse problems, some solvers in the literature were proposed by using numerical techniques and experimental testing. For instance, various publications in [1–10] recommended uses of the damping coefficient, stiffness, and external force for solving the inverse problems. Mode shape, frequency, displacement, and velocity at different times can also be used to estimate these properties successfully [11]. Huang [12] has employed a conjugate gradient method (CGM) to solve the nonlinear inverse vibration problem for the estimation of the time-dependent stiffness coefficient. Recently, Liu [13, 14] has developed a Lie-group shooting method to study the inverse vibration problem for estimating the time-dependent damping and stiffness coefficients and simultaneously derived a closed-form solution to estimate the parameters.

A complete review on the developments of useful methods for the realm of nonlinear system identification can be found in [14]. Reference [15] also proposed the idea of a force state mapping method which is a simple procedure and allows a direct identification of the restoring force for nonlinear mechanical systems. This idea was further extended in [16–18]. Recently, Namdeo and Manohar [19] have modified the force state mapping technique with two alternative schemes of functional representation: reproducing kernel particle method and kriging technique and estimating the parameters of nonlinear system from measured time histories of response under known excitations. Although this numerical method exhibits the capability to reproduce polynomials of specified order and has been applied to mechanical experiments, yet how to ensure numerical stability and avoid noise disturbance are not reported.

The aim of this paper is to develop a simple, multistep regularization algorithm with easy numerical implementation and versatility. A simple polynomial interpolation can be considered as a fit for the time history of displacement response under known excitations. Theoretically, if the order of the polynomial interpolation is large enough, the approximation displacement response can be as accurate as possible. However, in practical applications, it is not a good method to easily increase the order to obtain a highly accurate approximation. In fact, doing so will lead to a highly ill-posed matrix with a high-order function (the Vandermonde matrix), which has been described in [20]. To resolve the ill-posed linear problems for the Vandermonde system, Beckermann [21] and Li [22] claimed that the optimal condition number of the Vandermonde matrix could be expected. Nevertheless, the results of Beckermann [21] and Li [22] show that in the best possible cases, the condition numbers of the Vandermonde matrices still grow exponentially with the order of the interplant polynomial. Because of this, no one is interpolating the data by the high-order polynomials in the usual bases but rather in the Chebyshev polynomials. Hence, how to alleviate the ill-poseness owing to high-order function becomes one of the main tasks in this paper. First of all, we introduce the characteristic length (CL) of computational time into the high-order polynomial expansion to relieve the ill-conditioning of the resulting coefficient matrix of the polynomial expansion and then ensure numerical stability. This concept was first proposed to deal with the Laplace equation using a physical quantity [13, 14, 23, 24]. Recently, the CL has been successfully extended to deal with the Laplace equation and sloshing wave problems [25–27]. Although the CL can enhance the numerical accuracy for solving ill-posed linear matrix, it cannot avoid the effect of measured errors for parameters identification problems. Therefore, how to overcome the instability of the mathematical procedure is quite important. In addition, a small disturbance of measured data has to be considered in the numerical algorithm because they could cause an error identification of the parameter. In order to identify an accurate and stable solution for longer time scales, some special techniques have been used, including of the singular value decomposition (SVD), the SVD with a regularization parameter determined by the L-curve method, and sensitivity analysis. Despite these efforts, the stability problem remains unresolved. To thoroughly overcome these difficulties, this paper further adopts the CL combined with the natural regularization technique [28] to track ill-posed linear problems in numerical procedures. One advantage of this regularization method is that it can determine whether a solution exists for a linear system with the noisy level.

Apart from the current section, Section 2 describes the mathematical formulation of the characteristic time expansion method and introduces the numerical procedure of the matrix CGM. Section 3 gives several numerical examples, including of Duffing’s oscillator, Duffing’s oscillator with negative linear stiffness, Van der Pol’s oscillator, Bouc-Wen class model, and the seat model, to compare results of our method with the analytical solutions. Finally, some concrete conclusions are drawn in Section 4.

2. Basic Formulation

A second-order ordinary differential equation (ODE) for the equation of motion is expressed as

where represents the displacement of a system response; and are the external excitation and restoring force, respectively. In order to obtain , a trivial rearrangement of (1) gives

Here, can be obtained if the quantities, and , on the right-hand side of (2) are known. In general, it is easier to measure the displacement at some discrete sampling times than to directly measure velocities and accelerations. Therefore, if is denoted as the measured displacement, the differentiation of displacements can be expressed as follows:

This is however a set of index-three differential algebraic equations (DAEs) [29], which is hard to solve because the amplification of small errors and perturbations in the displacement cause severe numerically ill-posed condition.

2.1. The Characteristic Time Expansion Method

The polynomial interpolation is defined as the interpolation of a given set of data by a polynomial. In other words, given some data points, the aim is to find a polynomial which exactly goes through these points of data.

According to (3), the displacement can be expressed as a polynomial expansion:

where denotes each discrete time, denotes the displacement at each time, denotes the final time, and denotes the unknown coefficient. In many engineering applications, one wants to interpolate the data as accurate as possible. But this is limited by the interpolation of data with ()-order polynomials, where the resultant Vandermonde matrices are highly ill-conditioned as measured by the Lebesgue constant .

In this study, we introduce the characteristic length (CL) into (6) and express as follows

where denotes the CL. Differentiation of (7) yields velocity and acceleration and they are expressed as follows:

The polynomial expansion in (7) can be used to describe the displacement of a system. Hence, (7) can be expressed as a linear equation system with :

We denote the above equation by

where is the vector of unknown coefficients.

2.2. The Matrix CGM for Ill-Posed Linear System

When a matrix is ill-posed and measured data contains noisy disturbances, it is difficult to ensure the stability of the system using the conventional regularization techniques. Therefore, Liu et al. [28] proposed a natural regularization method, which proves that a solution exists when ill-posed matrix and noisy disturbances occur. This method can be described by the following matrix equation:

If is the inversion of , numerically, becomes a left-inversion of . Then we have

Let

Given , say , can be directly obtained because is a given matrix. Hence, we have

When (11) and (14) are combined together, they create an over-determined system to calculate . The over-determined system can be written as

where

is an matrix with . Multiplying (14) by yields an matrix equation:

Besides the primal system shown in (10), we need to solve the dual system with

Applying the operators in (17) to and utilizing the above equation, that is, , we can obtain

where .

Replacing the in (19) by , we have a similar equation for the primal system in (10)

where .

Finally, when of (20) is calculated by the CGM, the restoring force, velocity and acceleration can be obtained from (2) and (8).

3. Numerical Examples

Example 1. In this case, we consider a Duffing oscillator [29] and a second-order ODE to describe the forced vibration of a nonlinear structure by

where the parameters are fixed as , , and . The restoring force can be expressed as follows:

In order to identify the restoring force as a function of , a monotonic function of is required. In this example, is used to obtain the external force and is given by

To test the stability of the numerical method, the order of the polynomial and computational time are increased. The restoring force in the initial and final time changed very rapidly. To understand the CL effect, , , and are fixed. The maximum estimation error of with different CLs, shown in Figure 1, is smaller than .

Figure 1: For Example 1, showing the error of estimation with different CLs.

It can be seen that including the CL into this case is efficient to overcome an ill-posed matrix. Furthermore, by fixing , the exact solutions for velocity and acceleration can be determined. The numerical results are shown in Figures 2, 3, and 4. According to the numerical results, the maximum estimation errors of are found to be smaller than . We can find that applying the CL and matrix regularization method with the CGM can provide highly stable and accurate solutions. In order to further test the stability of the present method under different noise levels, we also consider

as an input into the estimation equations, where is a random number in , and is a noise level. With different noise levels %, 3%, 5%, and 10%, the computed profile of restoring forces is shown in Figure 5. Figure 5 also shows that the maximum estimated errors of are smaller than with noisy disturbances. We can find that the CL can effectively overcome numerical instability under the effect of the high-order function and large noise disturbances. Hence, we can see that the present method has a highly numerical accuracy and stability under the effect of the high-order function and large noise.

(a)

(b)

Figure 2: For Example 1, (a) comparing the exact solution and numerical one for velocity and (b) displaying the estimated error of velocity.

(a)

(b)

Figure 3: For Example 1, (a) comparing the exact solution and numerical one for acceleration and (b) displaying the estimated error of acceleration.

(a)

(b)

Figure 4: For Example 1, (a) comparing the exact solution and numerical one for restoring force and (b) displaying the estimated error of restoring force.

(a)

(b)

Figure 5: For Example 1, (a) comparing estimated and exact restoring forces under different noise level and (b) displaying the error of estimation with , 3%, 5%, and 10%, respectively.

Example 2. In this case, the Van der Pol oscillator is one of the nonlinear benchmark problem, and is given by

In this equation, is given by , and then, the external force can be obtained as

In this calculation, by fixing , , and , the numerical accuracy and stability of different parameters can be tested, including , and , , respectively. The maximum estimation error of shown in Figures 6 and 7 are smaller than.

(a)

(b)

Figure 6: For Example 2, (a) comparing exact solution and numerical one for restoring force with and under 4 seconds and (b) displaying the estimated error of restoring force.

(a)

(b)

Figure 7: For Example 2, (a) comparing exact solution and numerical one for restoring force with and under 4 seconds and (b) displaying the estimated error of restoring force.

To test the numerical stability of increasing the computational time by 10 seconds, the parameters are fixed as , and , , respectively. The maximum estimation errors of , shown in Figures 8 and 9, are smaller than . From the numerical solutions in Figures 6–9, it shows that the present method can keep the same numerical accuracy with the increase of the CL when the computational time increases.

(a)

(b)

Figure 8: For Example 2, (a) comparing exact solution and numerical one for restoring force with and under 10 seconds and (b) displaying the estimated error of restoring force.

(a)

(b)

Figure 9: For Example 2, (a) comparing exact solution and numerical one for restoring force with and under 10 seconds and (b) displaying the estimated error of restoring force.

This example demonstrates the results of fixing the parameters , , and under different noise levels with %, 3%, 5%, and 10%. The computed profile of is plotted in Figure 10. Figure 10(a) compares the restoring force with exact one, and the maximum estimation error of shown in Figure 10(b) is smaller than . From numerical result, we can find that the maximum error occurs because the displacement is equal to zero. That is, this present method can overcome the effect of the high-order function and large noise simultaneously. Therefore, it is found that the proposed method is accurate especially when noisy disturbances are encountered.

(a)

(b)

Figure 10: For Example 2, (a) comparing estimated and exact restoring forces and (b) displaying the error of estimation with , 3%, 5%, and 10%, respectively.

Example 3. The Bouc-Wen class model is one of the most widely used to efficiently describe smooth hysteretic behavior in engineering application. For a structural element described by the Bouc-Wen classical model, the restoring force is written as

where is a post- and preyield stiffness ratio, denotes natural frequency, and is an auxiliary variable that represents inelastic behavior. The evolution of is determined by an auxiliary ordinary differential equation, which can be written in the form of

where denotes the derivative of with respect to time and and are parameters that control the scale and sharpness of the hysteresis loops, respectively. Parameters, and , are used to control the shape of the hysteresis loop. In order to estimate the velocity and the restoring force of the Bouc-Wen model, the group preserving scheme (GPS), which was proposed by Liu [30], is adopted. The restoring force obtained by the GPS is referred to as the exact restoring force. We consider a system that has the parameter values as , , , , , , , , , and ; the initial condition of is given as (0.0, 0.0, 0.1). The exactly computed profile ofis plotted in Figure 11.

Figure 11: For Example 3, restoring force surface based on GPS.

To obtain using the present method, the parameters , , , and are fixed. The computed profile of at 40 to 42 seconds is plotted in Figure 12(a), and the maximum estimated error of , shown in Figure 12(b), is smaller than . Further, under a noise of %, the computed profile of at 40 to 42 seconds is plotted in Figure 13(a). The maximum estimated error of , shown in Figure 13(b), is smaller than . We see from Figures 12 and 13 that the maximum estimated error of still keep in the order of under a noise of%. That is, we can use the present method to achieve a more accurate and stable solution under a large noisy level.

(a)

(b)

Figure 12: For Example 3, (a) comparing estimated and exact restoring forces and (b) displaying the estimation error of restoring force.

(a)

(b)

Figure 13: For Example 3, (a) comparing estimated and exact restoring forces with and (b) displaying the estimation error of restoring force.

Example 4. As in Example 3, we consider the viscosity damping effect into the Bouc-Wen classical model and estimate the restoring force described as

where is the viscosity damping ratio. In this example, the parameter values are fixed as , , , , , , , , , , , and ; the initial condition of is given as . The exact computed profile of is plotted in Figure 14.

Figure 14: For Example 4, restoring force surface based on the GPS.

In this case, we use the same solver with the same parameters of Example 3. The computed profile of at 5 to 7 seconds is plotted in Figure 15(a). The maximum estimated error of , shown in Figure 15(b), is smaller than . Again, with a noise of %, the computed profile of at 5 to 7 seconds is plotted in Figure 16(a). The maximum estimated error of , shown in Figure 16(b), is smaller than . It can be seen in Figures 15 and 16 that the maximum errors are smaller than under a noise of %. Therefore, we conclude that for the smooth hysteretic behavior by the Bouc-Wen classical model, the accurate and stable solutions in Examples 3 and 4 are available when the proposed method is adopted.

(a)

(b)

Figure 15: For Example 4, (a) comparing estimated and exact restoring forces and (b) displaying the estimation error of restoring force.

(a)

(b)

Figure 16: For Example 4, (a) comparing estimated and exact restoring forces with and (b) displaying the estimation error of restoring force.

Example 5. In order to test the numerical stability of the CTEM used for the restoring force problem of discontinuous type, the vehicle seat problem is considered. When the seated human body is exposed to vertical vibration, the single-degree of freedom model can be used to describe its seat-person mathematical behavior, as shown in Figure 17, and is given by [31]

The parameters of the discontinuous typed vehicle seat model are given as , , , , , and. Here the external force is given by , and the parameters are given by , , , and , respectively. The computed profile of by and 201 is shown in Figure 18(a). The maximum estimated error of , shown in Figures 18(b) and 18(c), is smaller than . Numerical results show that this present method does not exhibit the numerical oscillation (Gibb’s phenomenon) when the high-order function is used. Hence, this example demonstrates that the present method has a high level of accuracy and stability for the restoring force problem of discontinuous type.

Figure 17: For Example 5, seat-person model of single-degree of freedom system [31].

(a)

(b)

(c)

Figure 18: For Example 5, (a) comparing estimated and exact restoring forces, (b) displaying the estimation error of restoring force with , and (c) displaying the estimation error of restoring force with .

4. Conclusions

In nonlinear mechanical system analysis, the inverse vibration problem is difficult to solve under the measured data with noise. This paper has successfully combined the CTEM with a natural regularization algorithm to determine the unknown restoring force. Due to inclusion of the CL to retain high accuracy and stability, the approximation method can avoid the numerical instability caused by a high-order polynomial function. In addition, when the measured data is contaminated by a large noise, the errors can be controlled by utilizing a natural regularization technique and increasing the CL. In summary, the presented method is an effective and convenient approach to solve the inverse vibration problems.

References

G. M. L. Gladwell, Inverse Problems in Vibration, vol. 9 of Monographs and Textbooks on Mechanics of Solids and Fluids: Mechanics. Dynamical Systems, Martinus Nijhoff Publishers, Dordrecht, The Netherlands, 1986.

G. M. L. Gladwell and M. Movahhedy, “Reconstruction of a mass-spring system from spectral data I: theory,” Inverse Problems in Engineering, vol. 13, pp. 179–189, 1997.

P. Lancaster and J. Maroulas, “Inverse eigenvalue problems for damped vibrating systems,” Journal of Mathematical Analysis and Applications, vol. 123, no. 1, pp. 238–261, 1987.

L. Starek and D. J. Inman, “On the inverse vibration problem with rigid-body modes,” Journal of Applied Mechanics, vol. 58, no. 4, pp. 1101–1104, 1991.

L. Starek and D. J. Inman, “A symmetric inverse vibration problem with overdamped modes,” Journal of Sound and Vibration, vol. 181, no. 5, pp. 893–903, 1995.

L. Starek and D. J. Inman, “A symmetric inverse vibration problem with overdamped modes,” Journal of Sound and Vibration, vol. 181, no. 5, pp. 893–903, 1995.

L. Starek and D. J. Inman, “A symmetric inverse vibration problem for nonproportional underdamped systems,” Journal of Applied Mechanics, vol. 64, no. 3, pp. 601–605, 1997.

S. Adhikari and J. Woodhouse, “Identification of damping: part 1, viscous damping,” Journal of Sound and Vibration, vol. 243, no. 1, pp. 43–61, 2001.

S. Adhikari and J. Woodhouse, “Identification of damping: part 2, non-viscous damping,” Journal of Sound and Vibration, vol. 243, no. 1, pp. 63–68, 2001.

M. Feldman, “Considering high harmonics for identification of non-linear systems by Hilbert transform,” Mechanical Systems and Signal Processing, vol. 21, no. 2, pp. 943–958, 2007.

G. Kerschen, K. Worden, A. F. Vakakis, and J.-C. Golinval, “Past, present and future of nonlinear system identification in structural dynamics,” Mechanical Systems and Signal Processing, vol. 20, no. 3, pp. 505–592, 2006.

C.-H. Huang, “A non-linear inverse vibratrion problem of estimating the time-dependent stiffness coefficients by conjugate gradient method,” International Journal for Numerical Methods in Engineering, vol. 50, no. 7, pp. 1545–1558, 2001.

C.-S. Liu, “Identifying time-dependent damping and stiffness functions by a simple and yet accurate method,” Journal of Sound and Vibration, vol. 318, no. 1-2, pp. 148–165, 2008.

C.-S. Liu, “A Lie-group shooting method for simultaneously estimating the time-dependent damping and stiffness coefficients,” Computer Modeling in Engineering and Sciences, vol. 27, no. 3, pp. 137–149, 2008.

S. F. Masri, A. G. Chassiakos, and T. K. Caughey, “Identification of nonlinear dynamic systems using neural networks,” Journal of Applied Mechanics, vol. 60, no. 1, pp. 123–133, 1993.

E. F. Crawley and A. C. Aubert, “Identification of nonlinear structural elements by force-state mapping,” AIAA Journal, vol. 24, no. 1, pp. 155–162, 1986.

E. F. Crawley and K. J. O'Donnell, “Identification of nonlinear system parameters in joints using the force-state mapping technique,” AIAA Paper 86–1013, 1986.

S. Duym, J. Schoukens, and P. Guillaume, “A local restoring force surface method,” International Journal of Analytical and Experimental Modal Analysis, vol. 11, pp. 116–132, 1996.

V. Namdeo and C. S. Manohar, “Force state maps using reproducing kernel particle method and kriging based functional representations,” Computer Modeling in Engineering and Sciences, vol. 32, no. 3, pp. 123–159, 2008.

I. Gohberg and V. Olshevsky, “The fast generalized Parker-Traub algorithm for inversion of vandermonde and related matrices,” Journal of Complexity, vol. 13, no. 2, pp. 208–234, 1997.

B. Beckermann, “The condition number of real Vandermonde, Krylov and positive definite Hankel matrices,” Numerische Mathematik, vol. 85, no. 4, pp. 553–577, 2000.

R.-C. Li, “Asymptotically optimal lower bounds for the condition number of a real vandermonde matrix,” SIAM Journal on Matrix Analysis and Applications, vol. 28, no. 3, pp. 829–844, 2006.

C.-S. Liu, “A meshless regularized integral equation method for Laplace equation in arbitrary interior or exterior plane domains,” Computer Modeling in Engineering and Sciences, vol. 19, no. 1, pp. 99–109, 2007.

C.-S. Liu, “A MRIEM for solving the laplace equation in the doubly-connected domain,” Computer Modeling in Engineering and Sciences, vol. 19, no. 2, pp. 145–161, 2007.

Y.-W. Chen, C.-S. Liu, and J.-R. Chang, “Applications of the modified Trefftz method for the Laplace equation,” Engineering Analysis with Boundary Elements, vol. 33, no. 2, pp. 137–146, 2009.

Y.-W. Chen, C.-S. Liu, C.-M. Chang, and J.-R. Chang, “Applications of the modified Trefftz method to the simulation of sloshing behaviours,” Engineering Analysis with Boundary Elements, vol. 34, no. 6, pp. 581–598, 2010.

Y.-W. Chen, W.-C. Yeih, C.-S. Liu, and J.-R. Chang, “Numerical simulation of the two-dimensional sloshing problem using a multi-scaling Trefftz method,” Engineering Analysis with Boundary Elements, vol. 36, no. 1, pp. 9–29, 2012.

C.-S. Liu, H.-K. Hong, and S. N. Atluri, “Novel algorithms based on the conjugate gradient method for inverting ill-conditioned matrices, and a new regularization method to solve ill-posed linear systems,” Computer Modeling in Engineering and Sciences, vol. 60, no. 3, pp. 279–308, 2010.

C.-S. Liu, “A Lie-group shooting method estimating nonlinear restoring forces in mechanical systems,” Computer Modeling in Engineering and Sciences, vol. 35, no. 2, pp. 157–180, 2008.

C.-S. Liu, “Cone of non-linear dynamical system and group preserving schemes,” International Journal of Non-Linear Mechanics, vol. 36, no. 7, pp. 1047–1068, 2001.

L. Wei and J. Griffin, “The prediction of seat transmissibility from measures of seat impedance,” Journal of Sound and Vibration, vol. 214, no. 1, pp. 121–137, 1998.