Energy Method to Obtain Approximate Solutions of Strongly Nonlinear Oscillators

Abstract

We introduce a nonlinearization procedure that replaces the system potential energy by an equivalent representation form that is used to derive analytical solutions of strongly nonlinear conservative oscillators. We illustrate the applicability of this method by finding the approximate solutions of two strongly nonlinear oscillators and show that this procedure provides solutions that follow well the numerical integration solutions of the corresponding equations of motion.

1. Introduction

This paper focuses on introducing a nonlinearization procedure to derive approximate solutions of strongly conservative nonlinear oscillators by replacing the system potential energy by equivalent representation forms that are based on cubic- or quartic-like polynomial expressions. The main motivation for finding an equivalent representation form of the potential energy terms is based on the different approaches that use equivalent restoring forces to derive approximate solutions of nonlinear oscillators. For instance, Yuste and Sánchez in [1, 2] used a cubication approach to replace the restoring forces by equivalent cubic polynomial expressions by applying the weighted mean-square and the harmonic balance methods. By using Chebyshev polynomial expansions, Beléndez and coworkers replaced the restoring forces by an equivalent cubic-like representation form that transformed the original equation of motion into a Duffing equation [3–5]. Then, Elas-Zúñiga and coworkers used this idea and considered additional Chebyshev polynomial expansion terms to transform the restoring forces into cubic-quintic polynomial representation forms that improved the accuracy of the approximate solutions of different types of strongly nonlinear oscillators [6, 7]. Here in this paper, instead of finding equivalent representation forms of the system restoring forces, we focus our attention on replacing the system potential energy by an equivalent energy expression that will be used to derive the corresponding approximate solutions. The accuracy of our proposed procedure will be illustrated by finding the approximate solutions of two nonlinear oscillators: (a) the cubic-quintic Duffing and (b) the finite extensibility nonlinear oscillator (FENO). We will show that our proposed approach provides solutions that follow well the numerical integration solution of the corresponding equations of motion. In fact, we have found that, for the case of the cubic-quintic oscillator, the nonlinear transformation form that is based on a quartic polynomial provides approximate angular frequency values whose maximum errors, when compared to the exact ones, do not exceed 0.004%. Similarly, we found for the FENO oscillator that when the initial oscillation amplitude 𝐴 approaches to one, the errors achieved do not exceed the percentage value of 0.3787%.

2. A Nonlinearization Energy Procedure

Conservative nonlinear oscillators are commonly modelled by equations of the form 𝑥̈+𝐹(𝑥)=0;𝑥(0)=𝑥₁₀,𝑥̇(0)=0, (1) where 𝐹(𝑥) is the conservative system restoring force and 𝑥₁₀ is the initial oscillation amplitude. Integration of (1) with respect to the displacement 𝑥 leads to its energy equation of motion 𝑥̇²+𝑉(𝑥)+𝑐₀=0. (2) Here, 𝑥̇ and 𝑉(𝑥) represent the velocity and the system potential energy terms, respectively, and 𝑐₀ is an integration constant whose value is determined from the given initial conditions. Now we assume that the potential energy term 𝑉(𝑥) can be replaced by an equivalent representation form which can be obtained by applying available techniques such as the nonlinearization approach [8], the Chebyshev polynomials of the first kind [3–5], the weighted mean-square method [2, 9], or other techniques available in the literature. In this paper, we adopt Chebyshev polynomial and used its general definition to write in equivalent form the potential energy terms [10] 𝑉(𝑥)=^𝑁_𝑛=0𝑏_𝑛𝑥₁₀𝑇_𝑛(𝑥)=𝑏₀2+^𝑁_𝑛=1𝑏_𝑛(𝑥)𝑇_𝑛(𝑥),−1≤𝑥≤1, (3) where 𝑏_𝑛=2𝜋⁺¹₋₁1√1−𝑥²𝑉(𝑥)𝑇_𝑛(𝑥)𝑑𝑥, (4) and the first four 𝑇_𝑛(𝑥) polynomials of the first kind are given by 𝑇₁(𝑥)=𝑥;𝑇₂(𝑥)=2𝑥²𝑇−1;₃(𝑥)=4𝑥³−3𝑥;𝑇₄(𝑥)=8𝑥⁴−8𝑥²+1. (5) Thus, the equivalent representation form of the system potential energy can be written as 𝑉(𝑥)≡𝑏₀2+𝑏₁𝑇₁(𝑥)+𝑏₂𝑇₂(𝑥)+𝑏₃𝑇₃(𝑥)+⋯+𝑏_𝑛𝑇_𝑛(𝑥)≈𝑐+𝑎₁𝑥+𝑎₂𝑥²+𝑎₃𝑥³+𝑎₄𝑥⁴+𝑎₅𝑥⁵+⋯+𝑎_𝑛𝑥^𝑛, (6) where the forms of the parameters 𝑐, 𝑎_𝑖, 𝑖=1,2,…,𝑛, depend on the Chebyshev coefficient terms. Therefore, the energy equation (2) can be equivalently represented by the following expression: 𝑥̇²+𝑎₀+𝑎₁𝑥+𝑎₂𝑥²+𝑎₃𝑥³+⋯+𝑎_𝑛𝑥^𝑛≈0, (7) where 𝑎₀≡𝑐₀+𝑐. However, this representation form (7) needs to be further transformed to a more tractable expression. To achieve this transformation, we can use the nonlinearization processes discussed in [6, 7]. Before we perform this transformation, we need to identify which polynomial order has to be used during the application of the nonlinearization processes, that is, a cubic-like or a quartic-like form. To better understand the implication of the choice of the polynomial order, we recall that the Duffing equation 𝑑²𝑥𝑑𝑡²+𝛼𝑥+𝛽𝑥³=0 (8) has an energy equation representation given as 𝑥̇²+𝛼𝑥²+𝛽2𝑥⁴+2𝑐₀=0, (9) where 𝛼, 𝛽, and 𝑐₀ are known parameter values. Therefore, it could be obvious that a quartic-like polynomial form has to be chosen to replace the potential energy terms of (7). However, we explore here two different options to replace the potential energy terms: (a) a cubic-like polynomial and (b) a quartic-like polynomial expression. If we use option (a), the nonlinearization processes will lead us to the following expression: 𝐹_𝑐(Δ,𝛿,𝜖)=^𝜎₀𝑎₀+𝑎₁𝑥+𝑎₂𝑥²+𝑎₃𝑥³+𝑎₄𝑥⁴+𝑎₅𝑥⁵+𝑎₆𝑥⁶+⋯+𝑎_𝑛𝑥^𝑛−𝛿𝑥−𝜖𝑥³−Δ²𝑑𝑥⟶min, (10) where the parameters 𝛿, 𝜖, and Δ can be found from 𝜕𝐹_𝑐(Δ,𝛿,𝜖)𝜕𝛿=0,𝜕𝐹_𝑐(Δ,𝛿,𝜖)𝜕𝜖=0,𝜕𝐹_𝑐(Δ,𝛿,𝜖)𝜕Δ=0. (11) Thus, the equivalent representation of (2) has the form 𝑥̇²+𝑉(𝑥)≡𝑥̇²+𝛿𝑥+𝜖𝑥³+Δ≈0. (12) However, if we adopt the quartic-like polynomial form, (2) can be transformed to 𝑥̇²+𝑉(𝑥)≡𝑥̇²+𝛿𝑥²+𝜖𝑥⁴2+2Δ≈0, (13) where 𝐹_𝑞(Δ,𝛿,𝜖)=^𝜎₀𝑎₀+𝑎₁𝑥+𝑎₂𝑥²+𝑎₃𝑥³+𝑎₄𝑥⁴+𝑎₅𝑥⁵+𝑎₆𝑥⁶+⋯+𝑎_𝑛𝑥^𝑛−𝛿𝑥²−12𝜖𝑥⁴−2Δ²𝑑𝑥⟶min, (14) with 𝜕𝐹_𝑞(Δ,𝛿,𝜖)𝜕𝛿=0,𝜕𝐹_𝑞(Δ,𝛿,𝜖)𝜕𝜖=0,𝜕𝐹_𝑞(Δ,𝛿,𝜖)𝜕Δ=0. (15) Of course, other possible polynomial representation forms could be assumed, however, we focus here on the usage of the cubic-like and the quartic-like polynomial forms to determine the angular frequency expressions of nonlinear oscillators under consideration. Notice that in (10) and (14), the value of 𝜎 must be appropriately chosen to ensure that the equivalent potential energy terms are qualitatively and quantitatively similar to those of the original equations of motion. In this sense, we have followed the approach developed by Cai and coworkers in [8].

We will next study the feasibility of our proposed nonlinearization energy approach to derive the equivalent energy equations of the cubic-quintic Duffing equation and then we will explore the solution of the finite extensibility nonlinear oscillator (FENO).

3. The Cubic-Quintic Duffing Oscillator

Although the exact solution of this nonlinear oscillator is well known [11], we choose its equation of motion to illustrate the applicability of our proposed mathematical technique to learn more about some possible issues that could arise during its implementation to solve nonlinear systems. This oscillator has an equation of motion given as 𝑦̈+𝛼𝑦+𝛽𝑦³+𝛾𝑦⁵=0,with𝑦(0)=𝑦₁₀,𝑦̇(0)=0, (16) where 𝑦 denotes the displacement of the system and 𝛼, 𝛽, and 𝛾 are system constant parameters. If we set 𝑥=𝑦/𝑦₁₀, then (16) can be written as 𝑥̈+𝐴𝑥+𝐵𝑥³+𝐺𝑥⁵=0,with𝑥(0)=1,𝑥̇(0)=0. (17) Here 𝐴=𝛼, 𝐵=𝛽𝑦²₁₀, and 𝐺=𝛾𝑦⁴₁₀. To apply our nonlinear energy approach, we first integrate (17) with respect to 𝑥. This step provides the following energy equation: 𝑥̇²+𝐴𝑥²+𝐵2𝑥⁴+𝐺3𝑥⁶+2𝑐₀=0, (18) where the integration constant 𝑐₀ has the value of 𝑐₀1=−2𝐵𝐴+2+𝐺3. (19) First, we follow option (a) and transform the potential energy terms of (18) into a cubic polynomial by following our proposed energy nonlinearization method. Then, we use (10) and (11) to write the potential energy in equivalent form as 𝑉(𝑥)=𝐴𝑥²+𝐵2𝑥⁴+𝐺3𝑥⁶+2𝑐₀≡𝛿𝑥+𝜖𝑥³+Δ, (20) where 𝛿=5𝐴𝜎−𝜎1256³7𝐵+9𝐺𝜎²,𝜖=35𝐴+154𝜎127𝐵𝜎+5𝐺𝜎³,Δ=2𝑐₀−𝐴𝜎²+27𝐵𝜎⁴+60𝐺𝜎⁶,42 (21) where 𝜎 is a parameter whose value satisfies (10). Thus, the equivalent energy transformation form of (18) is given as: 𝑥̇²+𝛿𝑥+𝜖𝑥³+Δ≈0. (22) To find the angular frequency of (22), we next proceed with its integration and find that 𝑡_𝑐=^𝑥_𝑥01−𝛿𝑥+𝜖𝑥³≡+Δ𝑑𝑥₁₀1−𝑥−𝑥₁𝑥−𝑥₂𝑥−𝑥₃𝜖𝑑𝑥, (23) where 𝑥₁=2(3)^1/3𝛿𝜖−2^1/39𝜖²Δ+3𝜖³4𝛿³+27𝜖Δ^22/36^2/3𝜖9𝜖²Δ+3𝜖³(4𝛿³+27𝜖Δ²)^1/3,𝑥₂=−2𝛿𝜖−3𝑖+√3+2^1/33^1/61+√×3𝑖9𝜖²Δ+3𝜖³4𝛿³+27𝜖Δ^22/3×22^2/33^5/6𝜖9𝜖²Δ+3𝜖³4𝛿³+27𝜖Δ^21/3−1,𝑥₃=−2𝛿𝜖3𝑖+√3+2^1/33^1/61−√×3𝑖9𝜖²Δ+3𝜖³4𝛿³+27𝜖Δ^22/3×22^2/33^5/6𝜖9𝜖²Δ+3𝜖³4𝛿³+27𝜖Δ^21/3−1, (24) with 𝑖=√−1. The integration of (23) provides the following time expression: 𝑡_𝑐=−2−1𝜖𝑥₁−𝑥₂×𝐹sin⁻¹𝑥₁−𝑥₂𝑥−𝑥₂,𝑥₃−𝑥₂𝑥₁−𝑥₂−𝐹sin⁻¹𝑥₁−𝑥₂𝑥₀−𝑥₂,𝑥₃−𝑥₂𝑥₁−𝑥₂. (25) Here 𝐹 represents the incomplete elliptic integral of the first kind with modulus (𝑥₃−𝑥₂)/(𝑥₁−𝑥₂). Thus, the approximate angular frequency value of the cubic-quintic oscillator is given as 𝜔_𝑛𝑐=2𝜋𝑇_𝑐, (26) where 𝑇_𝑐=4𝑡_𝑐 on 𝑥∈[0,1]. Before we assess the accuracy of our derived angular frequency expression (26), we next follow option (b) and write the potential energy as a function of a quartic polynomial expression 𝑉(𝑥)=𝐴𝑥²+𝐵2𝑥⁴+𝐺3𝑥⁶+2𝑐₀≡𝛿₁𝑥²+𝜖₁𝑥⁴2+2Δ₁. (27) Notice that the parameters 𝛿₁, 𝜖₁, and Δ₁ can be determined by using (14) and (15) which yields 𝛿₁=13333𝐴−5𝐺𝜎⁴₁,𝜖₁=11111𝐵+10𝐺𝜎²₁,Δ₁=1386𝑐₀+5𝐺𝜎⁶₁,1386 (28) where the value of 𝜎₁ must satisfy (14). Now the solution of 𝑥̇²+𝛿₁𝑥²+𝜖₁𝑥⁴2+2Δ₁≈0 (29) yields 𝑡_𝑞=^𝑥_𝑥01−𝛿₁𝑥²+𝜖₁𝑥⁴/2+2Δ₁≡𝑑𝑥_𝑥𝑥01−𝑥²−𝑥₁𝑥²−𝑥₂=𝜖/2𝑑𝑥−2𝐵𝑥₂𝐹sin⁻¹𝑥√𝑥₁,𝑥₁𝑥₂−𝐹sin⁻¹𝑥₀√𝑥₁,𝑥₁𝑥₂, (30) where 𝑥₁=−𝛿₁−𝛿²₁−4𝜖₁Δ₁𝜖₁,𝑥₂=−𝛿₁+𝛿²₁−4𝜖₁Δ₁𝜖₁. (31)(32) Thus, the approximate angular frequency expression is given as 𝜔_𝑛𝑞=2𝜋/𝑇_𝑞, where 𝑇_𝑞=4𝑡_𝑞. Also notice that (25) and (30) can be rewritten in different forms by putting 𝑥 as a function of time and by using the elliptic function identities. However, we will not elaborate on these solutions here.

We next compare our derived solutions with respect to the exact one of (17) found in [11]. Figure 1 illustrates the percentage error attained by comparing the exact angular frequency values of (17) with those obtained by using (26) plotted against the initial oscillation amplitude values. In this case, we use the parameter values of 𝐴=1, 𝐵=10, and 𝐺=1 with 𝜎=1.0854. We can see from Figure 1 that the maximum error values do not exceed 1.194%. We next use the quartic nonlinearization approximate angular frequency value given by (30) and compute the corresponding error values. One can notice from Figure 2, for the same system parameter values of 𝐴, 𝐵, and 𝐺, with 𝜎₁=1.13731, that now the maximum error values do not exceed 0.004%. This represents a great improvement when compared to the estimated values of 𝜔_𝑛𝑐. We believe that the improvement achieved by using (30) could be due to the usage of a fourth-order polynomial to approximate the potential energy expression which is closer to the order of the sextic polynomial that is used to derive, through the energy equation, the exact solution of (16). To further assess the precision of our proposed approach, we shall next explore the solution of a finite extensibility nonlinear oscillator (FENO).

Figure 1 

Percentage relative angular frequency error values curve versus the initial oscillation amplitude 𝑥10. The parameter values are 𝐴=1, 𝐵=10, 𝐺=1, and 𝜎=1.0854.

Figure 2 

Percentage relative angular frequency error values curve versus the initial oscillation amplitude 𝑥10. The parameter values are 𝐴=1, 𝐵=10, and 𝐺=1.

4. A Finite Extensibility Nonlinear Oscillator

This oscillator is used to describe chain dynamics of polymer-like materials subjected to high extensibility values [12]. Its equation of motion has the form 𝑑²𝑥𝑑𝑡²+𝑥1−𝑥²=0,𝑥(0)=𝐴,𝑥̇(0)=0, (33) where 𝐴 is the initial oscillation amplitude with 𝐴∈(0,1). If we introduce the transformation 𝑦=𝑥/𝐴, then (33) can be written as 𝑑²𝑦𝑑𝑡²+𝑦1−𝐴²𝑦²=0,𝑦(0)=1,𝑦̇(0)=0. (34) It is well known that when 𝐴→1, this oscillator tends to behave as an anharmonic one and thus, depending on the methods used, its approximate solutions to predict its dynamical response could diverge from the numerical one. See, for instance, [7–13]. The aim of this section focuses on exploring the applicability of our proposed nonlinearization energy approach to develop an approximate solution to (34) valid for the whole interval value of 𝐴.

First, we derive the energy expression of (34) by performing its integration with respect to 𝑦. This step leads to the following equation: 𝑦̇²−Log1−𝐴²𝑦²𝐴²+2𝑐₀=0, (35) where 𝑐₀=Log1−𝐴²2𝐴². (36) Since we need to write (34) in a polynomial-like form, then we follow our proposed approach and use the Chebyshev polynomial expansion (3) to get that 𝑦̇²+𝛼₁𝑦²+𝛽₁𝑦⁴+𝑐₁=0,𝑐₁=𝑐+𝑐₀, (37) where1𝑐=𝐴⁶8+3𝐴⁴−81−𝐴²+4𝐴²21−𝐴²−3−2𝐴⁴Log1/21−1−𝐴²,𝛼₁=−12𝐴⁴+8𝐴²9−51−𝐴²+641−𝐴²−1𝐴⁶,𝛽₁=88+𝐴⁴−81−𝐴²+4𝐴²1−𝐴²−2𝐴⁶. (38) We now adopt options (a) and (b) and determine the equivalent system potential energy. By using the cubic transformation form, we obtain that 𝑉𝑦=𝛼₁𝑦²+𝛽₁𝑦⁴+𝑐₁≡𝛿𝑦+𝜖𝑦³+Δ, (39) where 𝜎𝛿=−123𝛽₁𝜎²−5𝛼₁,7𝜖=54𝜎5𝛼₁+9𝛽₁𝜎²,1Δ=270270𝑐₁−10𝛼₁𝜎²+9𝛽₁𝜎⁴. (40) If we now follow option (b), the equivalent system potential energy is given as 𝑉𝑦=𝛼₁𝑦²+𝛽₁𝑦⁴+𝑐₁≡𝛿₁𝑦²+𝜖₁𝑦⁴2+2Δ₁, (41) for which 𝛿₁=𝛼₁, 𝜖₁=2𝛽₁, and Δ₁=𝑐₁/2.

Figure 3 illustrates the percentage error attained by comparing the angular frequency values of 𝜔_𝑛𝑐 and 𝜔_𝑛𝑞 with respect to the numerical integration ones plotted versus the initial oscillation amplitude 𝐴. It is interesting to note that when 𝐴 approaches to one, that is, 𝐴=0.9999999, the percentage error curves tend to the values of 0.3787% and 0.2040%, respectively. Figure 4 shows a comparison among our derived solutions and the angular frequency value obtained by Beléndez et al. [13] in which 𝜔_Bel=24−3𝐴²−5𝐴𝑐₁−30𝑐²₁𝐴,𝑐₁=8192𝐴³−12288𝐴⁵+6512𝐴⁷−1217𝐴⁹832768−65536𝐴²+52096𝐴⁴−19383𝐴⁶+2811𝐴⁸. (42)(43) It is clear from Figure 4 that when 𝐴 approaches to one, the highest percentage error is obtained from (42). Therefore, we can conclude that our derived angular frequency solutions provide error values that are remarkably lower than those previously reported in [7–13]. Furthermore, the maximum error values of both curves do not exceed 3.2%. Here the fitting parameter 𝜎 has the value of 1.01195.

Figure 3 

Percentage relative angular frequency error value curves versus the initial oscillation amplitude 𝐴. Here the fitting parameter 𝜎 has the value of 1.01195 and 0<𝐴≤0.9999999.

Figure 4 

Comparison of the percentage relative angular frequency error value curves versus the initial oscillation amplitude 𝐴.

We next focus on the qualitative behavior of our derived solution by plotting its amplitude-time response curve and compare this with its numerical integration solution. In this case, we have used (25) and (30) and considered half periods to ensure that our FENO system satisfies the initial conditions. In fact, we found that, for given values of 𝑦 bigger than one, the estimated values of 𝑡_𝑐 and 𝑡_𝑞 computed from (25) and (30), respectively, become imaginary. Therefore and in order to plot the amplitude-time response curves, the values of 𝑦 must satisfy the condition 0≤𝑦≤|1|. Figure 5 shows the amplitude-time response curves obtained from (25) and the numerical integration solution of (35). Notice that both solutions have the same qualitative behavior even though the assumed initial system oscillation amplitude 𝐴 has a value close to one; that is, 𝐴=0.999.

Figure 5 

Computed amplitude-time response curves of (34) and (25). Here the solid line represents the numerical integration solution of (34) while the black dashed line represents our derived energy approximate solution (25), by considering an initial oscillation amplitude value of 𝐴=0.999.

5. Conclusions

We have introduced a nonlinearization energy method to predict the qualitative and quantitative behavior of strongly nonlinear oscillators that is based on the nonlinear transformation of the system energy equation. We have assumed that the original system potential energy can be equivalently written as a cubic or quartic polynomial representation by using Chebyshev polynomials and the nonlinearization method. We have shown that this procedure leads to equations that can be exactly integrated by using elliptic functions. Then, we have used these equations to study the dynamical response of two nonlinear oscillators: one related to the cubic-quintic Duffing equation and the other related to a finite extensibility nonlinear oscillator (FENO). We have found that our proposed approach follows well the numerical integration solution of the corresponding equations of motion. In fact and for the case of the cubic-quintic oscillator, the nonlinear transformation approach based on a quartic polynomial provides approximate angular frequency values whose maximum errors do not exceed 0.004%. We believe that the improvement achieved by replacing the system potential energy of (18) by an equivalent potential energy form could be due to the usage of a fourth-order polynomial expression which is closer to the order of the sextic polynomial that is needed to derive, through the energy equation, the exact solution of (17).

Similarly, we found, for the FENO oscillator, that when the initial oscillation amplitude 𝐴 approaches to one, that is, when 𝐴=0.9999999, the errors achieved do not exceed the percentage value of 0.3787%. Furthermore, we have proved that the qualitative behavior of our derived solutions follows well that of the corresponding numerical integration solutions. Therefore, we can conclude that our nonlinearization energy approach provides results that describe well the dynamical response of these strongly nonlinear oscillators.

It is evident that our proposed approach could be used to study the dynamical response of other dynamical oscillators and that it can be expanded to obtain approximate solutions of dynamical system with two or more degrees of freedom.

Conflict of Interests

The authors declare that they have no conflict of interests with any mentioned entities in the paper.

Acknowledgments

This work was funded by Tecnológico de Monterrey, Campus Monterrey, through the Research Chair in Nanomaterials for Medical Devices and the Research Chair in Intelligent Machines. Additional support was provided by the European Union Seventh Framework Programme (FP7 PEOPLE-2009) under the Grant Agreement IRSES no. 247476 and from Consejo Nacional de Ciencia y Tecnología (Conacyt), Mexico.