# A New Method for Detecting the Time-Varying Nonlinear Damping in Nonlinear Oscillation Systems: Nonparametric Identification

Authors: T. S. Jang; Hyoungsu Baek; M. C. Kim; B. Y. Moon
Journal: Mathematical Problems in Engineering
 (2011)

---

## Abstract

This paper presents an original method that can be used for identifying time-varying nonlinear damping characteristics of a nonlinear oscillation system. The method developed involves the nonparametric identification, in which only the system responses, namely, displacement and velocity need to be known for the identification. However, the method is concerned with a Volterra-type integral equation of the first kind, which leads to an instability of numerical solutions. That is, the solutions identified lack stability properties. In order to overcome the difficulty, a stabilization technique is applied to the identification process. A numerical example comprising a highly nonlinear system is examined to demonstrate the workability of the proposed method for the time-varying damping identification.


---

## Body

## 1. Introduction

As a mathematical modeling, the identification of system damping has received increased attention in the fields of science and engineering over the last few decades. A proper modeling of nonlinear damping characteristics is important in many practical applications. For example, it is essential for the response prediction, real-time monitoring of structure, load identification, and so forth, in regard to nonlinear physical systems.Many identification methods have been investigated for identifying dampings. Iourtchenko et al. [1] successfully proposed an identification method for the case of parametric excitation. By using a measured stationary response, a procedure for the in-service identification of the damping characteristics was suggested on the basis of the stochastic averaging method of Iourtchenko and Dimentberg [2]. There are challenging issues associated with the modeling and identification of nonlinear damping systems [3–6].Even if significant progress has been made, many of challenging problems concerning the damping identification still await their solutions. Most of (nonlinear) damping identification is applied to autonomous systems: that is, the dampings we shall identify do not explicitly involve time. Recently, however, much attention has been paid to the time-varying damping identification of nonlinear dynamic models in various fields of engineering. For example, the estimation of time-varying damping is crucial in assessing the risk of negative aeroelastic damping in a wind turbine [7]. In this paper, we are concerned with the identification of such time-varying damping.We propose an original method for identifying time-varying nonlinear damping in nonlinear dynamical systems. The method proposed here is a novel one with remarkable improvements. First, based on the fact that the method is of nonparametric identification, any “priori” assumption of the functional form of time-varying damping is not required. Measured system responses of displacement and velocity are sufficient for the identification process of the method. Second, the method converts a measured displacement into a pseudodisplacement whose relationship with the time-varying damping is linear. Third, a numerical procedure of the method utilizes a stabilization technique in order to solve a problem of numerical instability occurred in the identification process. The stabilization technique is known as regularization, the idea of which is briefly to replace integral equations of the first kind by nearby second kind integral equations. Fourth, because the method does not depend upon small parameters at all, it overcomes the limitations and disadvantages of the technique of perturbation expansion. Finally, we demonstrate the workability of the method through a numerical experiment.


## 2. The Equation of Motion

A nonlinear oscillator with a time-varying nonlinear damping is considered, in which the governing equation is described as a second-order nonlinear differential equation:(2.1)mÿ+ky=λ(t)f(ẏ)+r(y),
where the symbols of m and k indicate constants. With αand β real values, an initial condition is imposed as(2.2)y(0)=α,ẏ(0)=β.
Equation (2.1) together with (2.2) constitutes an initial value problem for the motion of the oscillator. In this paper, we assume that the time-varying nonlinear damping of systems is characterized by a separable function, λ(t)f(ẏ), as in (2.1), where f(ẏ) denotes a general nonlinear function of the velocity ẏ and λ(t) its time-varying coefficient.We will transform the differential representation of (2.1) into an integral one, because sometimes mathematically an appropriate manipulation of an integral expression has advantages over that of differential one [8–14]. For example, integral operators appearing in integral expressions are in general stable, whereas differential operators of differential equations are usually unstable. Using a Green’s function method as well as a method of variation of parameters [15], we are able to find an integral relationship, equivalent to (2.1), as follows [9, 10]:(2.3)y(t)=αy1(t)+βy2(t)+∫0tK(t,τ)[λ(t)f(ẏ)+r(y)]dτ.
With Wronskian W=y1ẏ2-ẏ1y2, the kernel K in (2.3) is expressed as(2.4)K(t,τ)≡[y1(τ)y2(t)-y1(t)y2(τ)]mW(τ),
in which y1 and y2 are functions that satisfy the following equations:(2.5)mÿ1+ky1=0,y1(0)=1,ẏ1(0)=0,mÿ2+ky2=0,y2(0)=0,ẏ2(0)=1.


## 3. Integral Equation

Equation (2.3) is rewritten as(3.1)y(t)-αy1(t)-βy2(t)-∫0tK(t,τ)[r(y(τ))]dτ=∫0tK̃(t,τ)λ(τ)dτ,
where a new kernel K̃ is defined as(3.2)K̃(t,τ)≡K(t,τ)f[ẏ(τ)].
Since all the quantities on the left-hand side of (3.1) are known, (3.1) is an integral equation for λ(t). Let us call the left-hand side of (3.1) a pseudodisplacement η in this study, that is,(3.3)η(t)≡y(t)-αy1(t)-βy2(t)-∫0tK(t,τ)[r(y(τ))]dτ.
It then immediately follows that(3.4)η(t)=∫0tK̃(t,τ)λ(τ)dτ
from (3.1) and (3.3).The physical meaning of the integral equation (3.4) is that if we know the system’s responses of the displacement y and its velocity ẏ in advance, we can recover (or identify) the time-varying damping coefficient λ(t) in (2.1) by solving (3.4). In practice, conventional measurement devices enable the measuring of system’s responses to be performed in an accurate manner.Even though the oscillators with time-varying damping, considered herein, are “nonlinear” as exhibited in (2.1), the present problem to identify the oscillators’ time-varying damping turns out to be “linear”: that is, we establish a linear (integral) relationship between the pseudodisplacement η and the unknownλ(t) through (3.4).


## 4. Uniqueness

Prior to the discussion of a solution procedure for (3.4), physically it plays an important role to verify that the time-varying coefficient λ(t) in (3.4) can be identified (or determined) in a really unique manner. For the proof of uniqueness, we first remind that the integral equation of (3.4) is linear inλ(t). It thus suffices to show that the vanishing of the left-hand side of (3.4) (i.e., η≡0) implies λ≡0 in (3.4).Let us assume thatη in (3.3) vanishes. That is, (4.1)η(t)=y(t)-αy1(t)-βy2(t)-∫0tK(t,τ)[r(y(τ))]dτ=0,
or(4.2)y(t)=αy1(t)+βy2(t)+∫0tK(t,τ)[r(y(τ))]dτ.
Equation (4.2) is then a nonlinear integral equation that y should satisfy. Recalling that (2.1) is equivalent to (2.3), the differential equation which is equivalent to (4.2) is(4.3)mÿ+ky=r(y).
Comparing (2.1) and (4.3), we instantly arrive at(4.4)λ≡0.
This completes the proof. Thus, it has been shown that (3.4) has a unique solution. This means that we can recover the time-varying damping coefficient λ(t) in (2.1) in a unique way, once we know the measured responses of the displacement y and velocity ẏ.


## 5. Lack of Stability

Although we have found that the present (inverse) problem of determiningλ(t) is one-to-one, there still remains the question of the stability of the solution. That is, it is necessary to check whether the solution depends on the measured response data (displacement y and velocityẏ) in a continuous manner.According to the theory of first-kind integral equations, the “first” kind Volterra integral equation with a regular kernel such asK̃(t,τ) in (3.2) is ill-posed in the sense of stability [16–18]. Thus, the solution to (3.4) lacks stability properties. This means that the presence of noise in the measurement of y and ẏ can be greatly amplified in the numerical solution of λ(t), creating unreliable solutions.This (undesirable) phenomenon is in sharp contrast to the case of the usual well-posed integral equations of the “second” kind: conventional numerical methods are suitable for solving the usual well-posed equations, but they cannot be used to solve the ill-posed problem of (3.4) because of the lack of solution stability [16–18]. In fact, the direct (numerical) discretization of ill-posed systems such as (3.4) would result in ill-conditioned numerical systems with an extremely large condition number [8, 18–22].


## 6. Computingλ(t)

As mentioned in the previous section, even though many numerical methods have been developed for solving integral (or differential) equations, they are usually limited to well-posed (stable) problems. Thus, in order to solve the problems which lack stability properties, a special numerical procedure, suitable for the unstable problems, is required. To this end, we suggest the use of a stabilization technique, known as regularization. In this study, we use Tikhonov regularization to suppress the solution instability caused by the lack of stability properties.


### 6.1. Tikhonov Regularization

We briefly explain a procedure for Tikhonov regularization, the results from which will be presented at later section. According to Tikhonov [16], we construct a functional having an artificial dissipative function Ω for the solution for the integral equation of (3.4) for T>0,(6.1)F(λ)={∫0T|∫0tK̃(t,τ)λ(τ)dτ-η(t)|2dt}1/2+μ⋅Ω(λ),
where the constant μ is called the regularization parameter, which is a positive real number. In this study, the functional form of Ω in (6.1) is chosen as follows:(6.2)Ω(λ)=∫0T|d2λdt2|2dt.
We are able to minimize the (Tikhonov) functional F in (6.1) by using the orthogonal structure of an infinite dimensional Hilbert space. This leads to(6.3)μd4λdt4+L*L(λ)=L*(η).L in (6.4) denotes an integral operator, whose kernel is K̃(t,τ) in (3.2): that is,(6.4)L(ξ)=∫0tK̃(t,τ)ξ(τ)dτ,
for an arbitrary function ξ(t). The symbolic notations for L* in (6.3) is the adjoint operator of the L[8, 9]. Due to the role of the artificial function μ·Ω introduced in (6.1), we can arrive at the integrodifferential equation of (6.3), whereas the original equation remains the “first” kind integral equation of (3.4). This means that the solution for (6.3) is stable, because the integrodifferential equation (6.3) is known to be well-posed.


## 6.1. Tikhonov Regularization

We briefly explain a procedure for Tikhonov regularization, the results from which will be presented at later section. According to Tikhonov [16], we construct a functional having an artificial dissipative function Ω for the solution for the integral equation of (3.4) for T>0,(6.1)F(λ)={∫0T|∫0tK̃(t,τ)λ(τ)dτ-η(t)|2dt}1/2+μ⋅Ω(λ),
where the constant μ is called the regularization parameter, which is a positive real number. In this study, the functional form of Ω in (6.1) is chosen as follows:(6.2)Ω(λ)=∫0T|d2λdt2|2dt.
We are able to minimize the (Tikhonov) functional F in (6.1) by using the orthogonal structure of an infinite dimensional Hilbert space. This leads to(6.3)μd4λdt4+L*L(λ)=L*(η).L in (6.4) denotes an integral operator, whose kernel is K̃(t,τ) in (3.2): that is,(6.4)L(ξ)=∫0tK̃(t,τ)ξ(τ)dτ,
for an arbitrary function ξ(t). The symbolic notations for L* in (6.3) is the adjoint operator of the L[8, 9]. Due to the role of the artificial function μ·Ω introduced in (6.1), we can arrive at the integrodifferential equation of (6.3), whereas the original equation remains the “first” kind integral equation of (3.4). This means that the solution for (6.3) is stable, because the integrodifferential equation (6.3) is known to be well-posed.


## 7. Numerical Experiment

### 7.1. A Model Equation

In this section, we demonstrate the workability of the method proposed for the indirect measurement of time-varying damping. We begin with a model differential equation for the time-varying system of (2.1):(7.1)ÿ+λ(t)|ẏ|ẏ+y+0.1y3=0.
Recalling (3.2), we note that (7.1) is involved in the nonlinear function f in (3.2) which has the following form:(7.2)f[ẏ(τ)]=|ẏ(τ)|ẏ(τ).
Equation (7.2) is a typical functional form of the quadratic nonlinear damping arising in nonlinear dynamics of oscillations. Here, we list three cases of different λ(t)’s for the present numerical experiment as shown in Table 1. In case 1, a time-varying damping coefficient is set to λ(t)=0.1(1+e-0.25(t-7)2), which is bell-shaped centered at t=7. In contrast to such short transient period of the damping coefficient, case 2 has a periodic damping coefficient of λ(t)=0.1(1+0.5sin(t)). Case 3 has the combined feature of the cases 1 and 2, that is, a periodic oscillation within a slowly varying smooth envelope. The damping coefficients λ(t) of all the three cases are plotted in Figure 1.Table 1
Three test cases.
Caseλ(t)10.1(1+e-0.25(t-7)2)20.1(1+0.5sin(t))30.1(1+0.5sin(2πt))e-0.1(t-10)2Time-varying damping coefficients of three test cases.
(a)(b)(c)


### 7.2. System Responses

In order to solve the (inverse) problem in (3.4), measurements of displacement and velocity are needed as explained in Section 3. In general, such data is obtained by numerical solutions of the (direct) problem in (7.1). The system responses of the numerical solutions can be found by using various usual numerical integration schemes of nonlinear ordinary differential equations.In this study, the solution or measurement was obtained by using a second-order Runge-Kutta integration scheme, where the initial conditions are imposed as(7.3)y(0)=0,ẏ(0)=1.Small time steps and high-order Runge-Kutta integration methods are tested and we confirmed that the measurement, that is, numerical solution of the equation, is accurate in this study. In practice, the measurement of motion responses is always degraded to a certain extent by noise. Thus, it may be impossible to determine the left-hand side of (3.4) exactly. However, in this numerical experiment, we assume that the noise is of the order of discretization errors.As an example of a typical measurement, system responses are shown in Figures2 and 3. Due to positive system damping, the phase diagram shows a gradual decay of system energy (spiral-in in the phase plot). The final time T is set to 20 seconds. Although it has nothing to do with an inverse problem itself, numerically a large final time T with small time step will result in a large discretization numerical operator, say, a large matrix by a discretization of L in (6.3). In such situations, an iterative method like the Landweber regularization [9] can be very time-consuming and noniterative methods are preferred. This is a reason we present Tikhonov regularization in Section 6, which is a noniterative method.Time history of displacement and velocity (case 1).
(a)(b)Phase diagrams (case 1 in (a) and case 3 in (b)).
(a)(b)


### 7.3. Recoveringλ

λ(t) for three different regularization parameters of μ are illustrated in Figure 4 in which the Tikhonov regularization of (6.3) is used for the recovering. For the cases 2 and 3, the recovered λ(t)are depicted in Figures 5 and 6, respectively. Figures 4–6 show that the recovered results are in good agreement with the exact ones. But there is a small difference near the end points, which may be caused by nearly zero main diagonals involving a numerical discretization matrix for L in (6.3). When μ is too large, the smoothness is enforced so strongly that the recovered λ(t) cannot follow the sinusoidal (rapid) changes. When μ is too small, however, (3.4) is of a more dominating factor of minimization; hence, λ(t) follows the details of the damping coefficient. However, erroneous peaks appear at times when ẏ(t)is nearly zero due to the ill-posedness.Recovered damping coefficientλ(t) of case 1 (from (a) to (c) the regularization parameter).
(a)(b)(c)Recovered damping coefficientλ(t) of case 2 (from (a) to (c) the regularization parameter).
(a)(b)(c)Recovered damping coefficientλ(t) of case 3 (from (a) to (c) the regularization parameter).
(a)(b)(c)In this section, we demonstrate that the time varying damping can be covered when the regularization parameterμ is adjusted properly; in all cases, the time-varying coefficients are recovered successfully whether λ(t) has either short-period changes, periodic changes, or both (periodic changes with amplitude modulation).


## 7.1. A Model Equation

In this section, we demonstrate the workability of the method proposed for the indirect measurement of time-varying damping. We begin with a model differential equation for the time-varying system of (2.1):(7.1)ÿ+λ(t)|ẏ|ẏ+y+0.1y3=0.
Recalling (3.2), we note that (7.1) is involved in the nonlinear function f in (3.2) which has the following form:(7.2)f[ẏ(τ)]=|ẏ(τ)|ẏ(τ).
Equation (7.2) is a typical functional form of the quadratic nonlinear damping arising in nonlinear dynamics of oscillations. Here, we list three cases of different λ(t)’s for the present numerical experiment as shown in Table 1. In case 1, a time-varying damping coefficient is set to λ(t)=0.1(1+e-0.25(t-7)2), which is bell-shaped centered at t=7. In contrast to such short transient period of the damping coefficient, case 2 has a periodic damping coefficient of λ(t)=0.1(1+0.5sin(t)). Case 3 has the combined feature of the cases 1 and 2, that is, a periodic oscillation within a slowly varying smooth envelope. The damping coefficients λ(t) of all the three cases are plotted in Figure 1.Table 1
Three test cases.
Caseλ(t)10.1(1+e-0.25(t-7)2)20.1(1+0.5sin(t))30.1(1+0.5sin(2πt))e-0.1(t-10)2Time-varying damping coefficients of three test cases.
(a)(b)(c)


## 7.2. System Responses

In order to solve the (inverse) problem in (3.4), measurements of displacement and velocity are needed as explained in Section 3. In general, such data is obtained by numerical solutions of the (direct) problem in (7.1). The system responses of the numerical solutions can be found by using various usual numerical integration schemes of nonlinear ordinary differential equations.In this study, the solution or measurement was obtained by using a second-order Runge-Kutta integration scheme, where the initial conditions are imposed as(7.3)y(0)=0,ẏ(0)=1.Small time steps and high-order Runge-Kutta integration methods are tested and we confirmed that the measurement, that is, numerical solution of the equation, is accurate in this study. In practice, the measurement of motion responses is always degraded to a certain extent by noise. Thus, it may be impossible to determine the left-hand side of (3.4) exactly. However, in this numerical experiment, we assume that the noise is of the order of discretization errors.As an example of a typical measurement, system responses are shown in Figures2 and 3. Due to positive system damping, the phase diagram shows a gradual decay of system energy (spiral-in in the phase plot). The final time T is set to 20 seconds. Although it has nothing to do with an inverse problem itself, numerically a large final time T with small time step will result in a large discretization numerical operator, say, a large matrix by a discretization of L in (6.3). In such situations, an iterative method like the Landweber regularization [9] can be very time-consuming and noniterative methods are preferred. This is a reason we present Tikhonov regularization in Section 6, which is a noniterative method.Time history of displacement and velocity (case 1).
(a)(b)Phase diagrams (case 1 in (a) and case 3 in (b)).
(a)(b)


## 7.3. Recoveringλ

λ(t) for three different regularization parameters of μ are illustrated in Figure 4 in which the Tikhonov regularization of (6.3) is used for the recovering. For the cases 2 and 3, the recovered λ(t)are depicted in Figures 5 and 6, respectively. Figures 4–6 show that the recovered results are in good agreement with the exact ones. But there is a small difference near the end points, which may be caused by nearly zero main diagonals involving a numerical discretization matrix for L in (6.3). When μ is too large, the smoothness is enforced so strongly that the recovered λ(t) cannot follow the sinusoidal (rapid) changes. When μ is too small, however, (3.4) is of a more dominating factor of minimization; hence, λ(t) follows the details of the damping coefficient. However, erroneous peaks appear at times when ẏ(t)is nearly zero due to the ill-posedness.Recovered damping coefficientλ(t) of case 1 (from (a) to (c) the regularization parameter).
(a)(b)(c)Recovered damping coefficientλ(t) of case 2 (from (a) to (c) the regularization parameter).
(a)(b)(c)Recovered damping coefficientλ(t) of case 3 (from (a) to (c) the regularization parameter).
(a)(b)(c)In this section, we demonstrate that the time varying damping can be covered when the regularization parameterμ is adjusted properly; in all cases, the time-varying coefficients are recovered successfully whether λ(t) has either short-period changes, periodic changes, or both (periodic changes with amplitude modulation).


## 8. Discussion

We found that the formalism here has quite different characters from other related studies [8–14]; in the present study the kernel K̃(t,τ) in (3.4) depends on the data. More specifically, K̃(t,τ)in (3.2) includes system response data throughf[ẏ(τ)]=|ẏ(τ)|ẏ(τ); this may look naive and harmless. However, the velocity ẏ in (3.2) may be close to zero, which may induce the unwelcome small difference near the end points in Figures 4–6, just mentioned above. Moreover, this velocity data is contaminated with a noise. Such data dependency of kernel causes an additional difficulty in the regularization, giving rise to sharp peaks in the recoveredλ(t).


## 9. Conclusion

The basic idea underlying the present study can be summarized as a fundamental question. Is it mathematically possible to identify (or recover) the time-varying damping characteristic of a nonlinear oscillation system, if we obtain “two” sets of motion response data of measured displacement and velocity? As a first step, we have answered the above question by assuming that the time-varying nonlinear damping of the system is separable: that is, it is expressed as a multiplication of a time-varying coefficientλ(t) by a nonlinear function f(ẏ). Based on the separability assumption, we present a general method for identifying the time-varying coefficient λ(t) of the system. Although a highly nonlinear system is provided as an example, the proposed method achieves a reasonably good approximation for λ in a stable manner.



---

*Source: 749309-2011-11-16.xml*