An Analytical Approximation Method for Strongly Nonlinear Oscillators

Abstract

An analytical method is proposed to get the amplitude-frequency and the phase-frequency characteristics of free/forced oscillators with nonlinear restoring force. The nonlinear restoring force is expressed as a spring with varying stiffness that depends on the vibration amplitude. That is, for stationary vibration, the restoring force linearly depends on the displacement, but the stiffness of the spring varies with the vibration amplitude for nonstationary oscillations. The varied stiffness is constructed by means of the first and second averaged derivatives of the restoring force with respect to the displacement. Then, this stiffness gives the amplitude frequency and the phase frequency characteristics of the oscillator. Various examples show that this method can be applied extensively to oscillators with nonlinear restoring force, and that the solving process is extremely simple.

1. Introduction

Many methods have been developed for the study of strongly nonlinear oscillators and several of them have been used to find approximate solutions to nonlinear oscillators [1–11]. Due to the diversity of nonlinearity, some methods are limited by the parameters of the equation or invalid for the forced oscillator even for an oscillator with nonlinear restoring force. Because of this, most of the methods directly seek the harmonic component while neglecting the character of the nonlinearity. In fact, nonlinear restoring forces can be classified into only two kinds: softening and hardening springs. For the former/latter, the free-vibration frequencies decrease/increase along with the amplitude. The varying rate of the free-vibration frequency is decided by the restoring force that is characterized by its function and derivatives. The present approach proposes a method that directly relates the frequency and amplitude to the character of the restoring force. The preliminary exploration of this method solved the amplitude-frequency characteristics of hardening and piecewise linear spring [11]. Here, some mathematical explanation on the reasonability and generality of this method are introduced.

2. Method

It is well known that the frequency of the free oscillator with nonlinear restoring force is dependent on the amplitude. And the ways of frequency depending on the amplitude are different for hardening and softening springs. To quantitatively describe the characteristics of hardening and softening springs, we introduce the first and second averaged derivatives of restoring force 𝑓(𝑥) with respect to the displacement as follows:𝑆(𝐴)=[𝑓(𝐴)−𝑓(0)]𝐶(𝐴)=𝑓(𝐴)−𝑓(𝐴/2)𝐴,or𝐶(𝐴)𝐴=/24[𝑓(𝐴)−−𝑓2𝑓(𝐴(𝐴/2𝐴/)2/−2)+𝑓𝑓(0()0)]/(𝐴/2),𝐴²,(2.1) where 𝐴 is the vibration amplitude. The first averaged derivative, 𝑆(𝐴), describes the maximum restoring force for given vibration amplitude; and the second derivative, 𝐶(𝐴), indicates the nonlinearity of the spring. These two averaged derivatives give the characteristics of restoring force roughly, but directly and simply. The second derivative is zero in the case of linear spring, and it takes positive/negative for hardening/softening spring. From the view point of geometry, the first and second derivatives are the slope of line OP and slope varying rate of the broken line OCP, respectively, as shown in Figure 1.

Figure 1 

The geometric characterizing of the nonlinear restoring force.

From Figure 1, it can be seen that the triangle OCP transforms along with point 𝑃, namely, the vibration amplitude 𝐴. Based on the triangle, a line OL is constructed with the slope 𝑘(𝐴) defined by the following equation:𝑘(𝐴o)r𝐴𝑘/2(−𝐴)𝑓=(𝐴2/[2𝑓)=(𝐴𝑓)+(𝐴𝑓)(−𝐴𝑘/2()𝐴])𝐴,3𝐴.(2.2) Equation (2.2) shows that 𝑘(𝐴) is dependent on the amplitude 𝐴, and difference between 𝑓(𝑥) and 𝑘(𝐴)𝑥 is the same at 𝑥=𝐴/2 (point C) and 𝑥=𝐴 (point 𝑃), as shown in Figure 1. Here, 𝑘(𝐴)𝑥 can be taken as a linear restoring force for given amplitude, or a linear restoring force with varied stiffness 𝑘(𝐴). This does not exactly approximate the real restoring force, but gives the dependence of frequency on the amplitude simply. After replacing the nonlinear restoring force with 𝑘(𝐴)𝑥, the angular frequency of the nonlinear oscillator can be calculated as𝜔(𝐴)=√𝑘(𝐴).(2.3)

Now, the validity of (2.3) is preliminarily examined with oscillators that have single-term positive-power nonlinearity written as follows:𝑥̈+sign(𝑥)|𝑥|^𝑏=0.(2.4) The initial conditions are taken as𝑥(0)=𝐴,𝑥̇(0)=0.(2.5) Equation (2.4) is a conservative nonlinear oscillatory system having a rational form for the nondimensional restoring force. It has been demonstrated that all the curves in the phase-space corresponding equation (2.4) are closed, and all motions for arbitrary initial conditions give periodic solutions [2]. Considering (2.2)–(2.5), we have𝜔(𝐴)=√𝑘(𝐴)=231+(0.5)^𝑏𝐴^(𝑏−1)/2.(2.6) For 𝑏=1, (2.6) is degenerated to the linear case, 𝜔=1, and for 𝑏=3 it gives𝜔(𝐴)=34𝐴.(2.7) This is the same result as the first approximation of harmonic balance method (HBM) [3]. As an example, for 𝑏=11/4, (2.6) simply gives𝜔(𝐴)=√𝑘(𝐴)=231+21^11/4𝐴^7/8=0.8751𝐴^7/8.(2.8) To compare with the numerical result, we take 𝐴=5, the angular frequency calculated by (2.8) is 3.5781, and the numerical result is 3.5279 with the ODE solver Ode45 of Matlab. The error is about 1.37%. The “Harmonic Balance Fourier” method can give the exact result, but in terms of Gamma functions [7].

Numerical simulation shows that the frequency error calculated with (2.6) is almost linearly increased along with the exponent; this is the cost for simplicity. Here, the error is corrected by the following correction coefficient 𝛼 that is determined by a least squares fit with simulation:𝛼=⎧⎪⎨⎪⎩110.,.0946−+00..0024𝑏𝑏,,𝑏𝑏1><≤21𝑏,.≤2,(2.9) In the case of the multiterm powers, some powers may have negative coefficients, the correction coefficient is changed to𝛼=⎧⎪⎨⎪⎩011.,.9064+−00..0024𝑏𝑏,,𝑏𝑏1><≤21𝑏,.≤2,(2.10) Then, the angular frequency equation (2.6) is corrected to𝜔(𝐴)=𝛼√𝐾(𝐴),(2.11) and the corrected stiffness is𝑘̃(𝐴)=𝛼²𝑘(𝐴)=2𝛼²31+12^𝑏𝐴^𝑏−1.(2.12) For the restoring force consisting of several power terms, each term is dealt with individually to get its stiffness 𝑘̃_𝑏(𝐴), and then the total stiffness is to take their sum because 𝑘̃_𝑏(𝐴)𝑥 is linear in 𝑥. For the restoring force in the forms other than power terms, if it can be approximated or estimated by powers, for example 𝑎𝑥^𝑏, then we directly use (2.2) to get its stiffness and correct it according to the corresponding exponent b; if it can be expanded in power series, then it is dealt with as the multi-term powers.

3. Examples

The method presented here is different to the conventional ones, the main advantages are no limitation on the nonlinear forms of restoring force and the simple calculation. As an example, we constructed the following oscillator only for examining the method:𝑥̈+𝑥+2sign(𝑥)|𝑥|^7/4+3sign(𝑥)|𝑥|^19/6=0.(3.1) Since the fractional exponents, it is difficult or complex for the conventional method to get the dependency of frequency on amplitude. Using (2.12), we can directly write out the stiffness as𝑘̃_7/4(𝐴)=43𝛼^27/41+(0.5)^7/4𝐴^3/4,𝑘̃_19/6(𝐴)=2𝛼^219/61+(0.5)^19/6𝐴^13/6.(3.2) The total stiffness of the restoring force and the angular frequency are respectively,𝑘̃(𝐴)=1+𝑘̃_7/4(𝐴)+𝑘̃_19/6(𝐴),𝜔(𝐴)=𝑘̃(𝐴).(3.3) From (2.10), the correction coefficients are 𝛼_7/4=1 and 𝛼_19/6=0.9767. For 𝐴=5, 2, and 1, (3.2) gives, respectively𝜔(5)=8.7233,𝜔(2)=3.6645,𝜔(1)=2.2023.(3.4) The numerical results are 𝜔(5)=8.7013,𝜔(2)=3.6679,and𝜔(1)=2.2108 solved by solver Ode113 of Matlab. Even if this equation can be solved by other methods, the procedure is extremely complex. The main trouble is caused by the fractional exponent of the power.

The second example is the antisymmetric constant force oscillator [6, 7]:𝑥̈+sign(𝑥)=0.(3.5) It is equivalent to 𝑏=0 in (2.4). Using (2.12) and 𝛼=0.96 calculated from (2.9), we get𝑘̃(𝐴)=𝛼²𝑘(𝐴)=2×0.96²31+(0.5)⁰𝐴⁰⁻¹,𝜔(𝐴)=𝑘̃(𝐴)=1.1085√𝐴.(3.6) The exact result is 1.1107/√𝐴.

The third example is following oscillator [5]:𝑥̈+𝑥³1+𝑥²=0.(3.7) For rough approximation, using (2.2), we have𝑘(𝐴)=2𝐴²31/1+𝐴²+1/8+2𝐴²,𝜔(𝐴)=√𝑘(𝐴).(3.8) Although the restoring force is not in power form, it approaches 𝑥³ and 𝑥 for |𝑥|≪1 and |𝑥|≫1, respectively. Therefore, in the procedure of correcting the stiffness or frequency in (3.8), we take the nonlinear function in (2.10) as 𝑥³ and 𝑥 for |𝑥|<1 and |𝑥|≥1, respectively. In other words, we correct the frequency in (3.8) according to 𝑏=3 and 𝑏=1 for |𝐴|<1 and |𝐴|≥1, then we have 𝛼=0.98 and 𝛼=1, respectively. The comparison to the exact result 𝜔_𝑒 is listed in Table 1.

Table 1 

The exact and approximated frequency.

The last example in this section is the pendulum:𝑥̈+sin𝑥=0.(3.9) For 𝑥<0.5𝜋, the restoring force can be directly approximated with the stiffness𝑘(𝐴)=23𝐴sin𝐴+sin𝐴2.(3.10) The error of frequency is less than 0.45%, and it is increased to 2.67% for 𝑥=0.7𝜋. To get accurate result for large amplitude, we expand the restoring force assin𝑥≅𝑥−𝑥³6+1𝑥2⁵0−50𝑥4⁷0.(3.11) For nonlinear terms in the right hand, we write out their stiffness by using (2.12):𝑘̃₃(𝐴)=−0.125𝛼²³𝐴²,𝑘̃₅(𝐴)=0.0057𝛼²⁵𝐴⁴,𝑘̃₇(𝐴)=−1.333×10⁻⁴𝛼²⁷𝐴⁶.(3.12) The correction coefficient 𝛼₅ is 0.94 calculated by (2.9), while 𝛼₃ and 𝛼₇ are 1.02 and 1.1 calculated with (2.10), respectively. Then the total stiffness and angular frequency are obtained as𝑘̃(𝐴)=1+𝑘̃₃(𝐴)+𝑘̃₅(𝐴)+𝑘̃₇(𝐴).(3.13) The angular frequency depending on the amplitude is𝜔(𝐴)=√1−0.1301𝐴²+0.0051𝐴⁴−0.0001613𝐴⁶.(3.14) The amplitude-frequency diagram of (3.14) is plotted in Figure 2.

Figure 2 

The diagram of amplitude 𝐴 depending on frequency 𝜔 for pendulum.

For comparing with the numerical results we use respectively, the ODE solver Ode45, Ode113, and numerical integrating the following integration from  (3.9)𝜋2𝜔=^𝐴0𝑑𝑥√2(cos𝑥−cos𝐴).(3.15) Since the solver Ode45 and Ode113 do not give the same results, and Ode113 matches the result of numerical integration, therefore we take and plot the result of the latter in Figure 2, as denoted by asterisks *.

4. Examples of Forced Vibration

From the examples in the last section, it can be seen that the fundamental component of the free vibration can be given by the equivalent stiffness 𝑘(𝐴). Now, the validity to approximate the nonlinear restoring forces with 𝑘(𝐴)𝑥 is examined for the forced oscillator. The following equations are forced oscillators with hardening/softening spring, respectively.𝑥̈+2𝑐𝑥̇+𝑥±13sign(𝑥)|𝑥|^2.5=𝐹sinΩ𝑡.(4.1) In (4.1) the positive/negative signs of the nonlinear term correspond to hardening/softening spring, respectively. They are dealt with in the same way except the calculation of the correction coefficient. At first, suppose that (4.1) has the following stationery response:𝑥=𝐴sin(Ω𝑡+𝛽).(4.2) Then according to the amplitude 𝐴, the nonlinear part of the restoring force in (4.1) can be approximated as𝑘̃_2.5(𝐴)=±21+0.5^2.59𝛼^22.5𝐴^1.5.(4.3) Then adding the linear restoring force, we get the total stiffness:𝑘̃(𝐴)=1±21+0.5^2.59𝛼^22.5𝐴^1.5.(4.4) Then (4.1) can be rewritten as𝑥̈+2𝑐𝑥̇+𝑘̃(𝐴)𝑥=𝐹sinΩ𝑡.(4.5) Equation (4.5) formally is a linear equation for stationery motion, the constructed stiffness 𝑘̃(𝐴) keeps (4.5) the same fundamental harmonic response as (4.1). Substitution of (4.2) into (4.5) yields−𝐴Ω²sin(Ω𝑡+𝛽)+2𝑐𝐴Ωcos(Ω𝑡+𝛽)+𝑘̃(𝐴)𝐴sinΩ𝑡=𝐹sinΩ𝑡.(4.6) From (4.6), the amplitude and phase angle depending on the forcing frequency are obtained as𝐴²Ω²−tan𝑘̃𝛽(𝐴=)²+2(𝑐2Ω𝑐𝐴Ω)²=𝐹²,Ω²−𝑘̃(𝐴).(4.7) For the hardening and softening springs, the amplitude-frequency equations are respectively,𝐴𝐴²₂Ω²−1−0.25889𝐴^1.5²+(2𝑐𝐴Ω)²=𝐹²,Ω²−1+0.26412𝐴^1.5²+(2𝑐𝐴Ω)²=𝐹².(4.8) The amplitude frequency diagrams for hardening springs with 𝑐=0.07 and 𝐹=1.15, and for softening one with 𝑐=0.07 and 𝐹=0.15 are plotted in Figure 3.

(a)

(b)

Figure 3 

The amplitude depending on frequency 𝜔 for hardening (a) and softening (b) springs.

Equation (4.8) gives both the stable and unstable periodic solutions; here we do not discuss them further.

5. Conclusion

For oscillator nonlinear restoring force, the stiffness depending on amplitude is constructed with the restoring forces at displacements 0, 𝐴/2 and 𝐴; then the stiffness is used to evaluate the frequency for free vibration and to find the response for forced vibration. The proposed method is geometrically explained by means of the first and second averaged derivatives of the restoring force. Various examples show that the solution can be approximated simply and accurately, no matter how complex the restoring force function is.

Acknowledgment

This paper was partially supported by the National Natural Science Foundation of China under Grant no. 11172018.