Solution of Excited Non-Linear Oscillators under Damping Effects Using the Modified Differential Transform Method

Abstract

The modified differential transform method (MDTM), Laplace transform and Padé approximants are used to investigate a semi-analytic form of solutions of nonlinear oscillators in a large time domain. Forced Duffing and forced van der Pol oscillators under damping effect are studied to investigate semi-analytic forms of solutions. Moreover, solutions of the suggested nonlinear oscillators are obtained using the fourth-order Runge-Kutta numerical solution method. A comparison of the result by the numerical Runge-Kutta fourth-order accuracy method is compared with the result by the MDTM and plotted in a long time domain.

1. Introduction

Our concern in this work is to give semi-analytic solutions of excited Duffing and excited van der Pol oscillators under damping effect which are given in the forms

$$\begin{array}{cc}\hfill \text{Duffing}\text{equation}& \frac{d^{2}x}{d{t}^2}+\eta \frac{dx}{dt}+{\omega }^{2}x+\alpha x^3=Asin(\Omega t)\hfill \end{array}$$

(1a)

$$\begin{array}{cc}\hfill \text{van}\text{der}\text{Pol}\text{equation}& \frac{d^{2}x}{d{t}^2}-ϵ(1-x^2)\frac{dx}{dt}+{\omega }^{2}x=Asin(\Omega t)\hfill \end{array}$$

(1b)

where x is the position coordinate which is a function of the time t, ω is the system’s natural frequency, η is a scalar parameter indicating the damping factor in Duffing equation, ϵ the nonlinearity and strength of the damping in van der Pol equation, respectively. α is a nonlinear parameter factor, A and Ω are the forcing amplitude and frequency, respectively.

Those two considered nonlinear oscillators have received remarkable attention in recent decades due to the variety of their engineering applications. For example, Duffing Equation (1a) used in studying the magneto-elastic mechanical systems [1], nonlinear vibrations of beams and plates [2,3] and vibrations induced by fluid flow [4] which are modeled by the nonlinear Duffing equation.

On the other hand, during the first half of the twentieth century, Balthazar van der Pol pioneered the fields of radio and telecommunications [5,6,7,8,9,10]. In an era when these areas were much less advanced than they are today, vacuum tubes were used to control the flow of electricity of transmitters and receivers. Simultaneously with Lorenz, Thompson, and Appleton, van der Pol experimented with oscillations in a vacuum tube triode circuit and concluded that all initial conditions converged to the same orbit of a finite amplitude. Since this behavior is different from the behavior of solutions of linear equations, van der Pol proposed a nonlinear differential Equation (1b) without excitation force, i.e., $A=0$, commonly referred to as the (unforced) van der Pol equation [8], as a model for the behavior observed in the experiment. In studying the case $\eta <<1$, van der Pol discovered the importance of what has become known as relaxation oscillations [8].

The most common methods for constructing approximate analytical solutions to the nonlinear oscillator’s equations are the perturbation methods [11]. These methods include the harmonic balance method, the elliptic Lindstedt-Poincaré method [11,12,13]. The Krylov Bogoliubov Mitropolsky method [14,15], the averaging [11,12,13,14,15,16] and multiple scales method [12] are widely used to obtain approximate solutions of nonlinear oscillators. A general common factor to all of these methods is that they solve weakly the nonlinear systems by using perturbation techniques to reduce the system into simpler equations which transform the physical problem into a purely mathematical one, for which a solution is readily available.

This work is the derivation to obtain approximate analytical oscillatory solutions for the nonlinear oscillator Equations (1a) and (1b) with initial conditions $x\left(0\right)=a$ and $\dot{x}\left(0\right)=b$ using the modified differential transform method. This is a powerful method for solving linear and nonlinear differential equations. This method was at first used as differential transform method in the engineering domain by Zhou [17] and in fluid flow problems [18,19,20,21,22,23,24,25]. The differential transform solution diverges by using finite number of terms. To solve this problem the modified differential transform method [26,27,28,29,30] was developed by combining the differential transform method (DTM) with the Laplace transform and Padé approximant [31] which can successfully predict the solution of differential equations with finite numbers of terms [32,33].

2. Differential Transform Method

A brief explanation of the differential transform method (DTM) is given in [27,28,29,30,31]; for an analytic function $x\left(t\right)$ in domain G, which can be represented by a power series around any arbitrary point in this domain. The differential transform of $x\left(t\right)$ is defined as follows:

$$X\left(k\right)=\frac{1}{k!}{\left[\frac{d^{k}x\left(t\right)}{d{t}^k}\right]}_{t=0}$$

(2)

In Equation (2), $x\left(t\right)$ is the original function and $X\left(k\right)$ is the transformed corresponding function. The inverse transform of $X\left(k\right)$ is defined as

$$x\left(t\right)=\sum _{k=0}^{\infty }X\left(k\right)t^k$$

(3)

Combining Equations (2) and (3), we obtain the following equation

$$x\left(t\right)=\sum _{k=0}^{\infty }\frac{t^k}{k!}{\left[\frac{d^{k}x\left(t\right)}{d{t}^k}\right]}_{t=0}$$

(4)

In applications, $x\left(t\right)$ takes finite number of terms and Equation (4) can be written as

$$x\left(t\right)=\sum _{k=0}^N\frac{t^k}{k!}{\left[\frac{d^{k}x\left(t\right)}{d{t}^k}\right]}_{t=0}$$

(5)

Some basic transformation rules of the differential transform method which are used in this work are tabulated in

Table 1. The fundamental operations of the differential transform method (DTM).

Original Function	Transformed Function
$\alpha u\left(t\right)\pm \beta v\left(t\right)$	$\alpha U\left(k\right)\pm \beta V\left(t\right)$
$u\left(t\right)v\left(t\right)$	${\sum }_{t=0}^{k}U\left(l\right)V(k-l)$
$u\left(t\right)v\left(t\right)w\left(t\right)$	${\sum }_{s=0}^k{\sum }_{m=0}^{k-s}U\left(s\right)V\left(m\right)W(k-s-m)$
${\displaystyle \frac{d^{m}u\left(t\right)}{d{t}^m}}$	${\displaystyle \frac{(k+m)!}{k!}U(k+m)}$
$exp\left(t\right)$	${\displaystyle \frac{1}{k!}}$
$sin(\omega t+\alpha )$	${\displaystyle \frac{{\omega }^k}{k!}sin(k\pi /2+\alpha )}$
$cos(\omega t+\alpha )$	${\displaystyle \frac{{\omega }^k}{k!}cos(k\pi /2+\alpha )}$

.

Now, we will apply the MDTM to obtain semi-analytic solutions for forced nonlinear Duffing and van der Pol differential equations under damping effect.

3. Forced Duffing Oscillator under Damping Effect

Consider a nonlinear differential Equation (1a) which describes the forced Duffing oscillator with damping effect and initial conditions $x\left(0\right)=a$ and $\dot{x}\left(0\right)=b$. The differential transform of this equation gives the recurrence relation

$$\begin{array}{c}\hfill (k+2)(k+1)X(k+2)+{\omega }^{2}X\left(k\right)+\eta (k+1)X(k+1)\\ \hfill +\alpha \sum _{n=0}^k\sum _{m=0}^{n}X\left(m\right)X(n-m)X(k-n)-\frac{A{\Omega }^k}{k!}sink\pi /2=0\end{array}$$

(6)

$$\begin{array}{c}\hfill X\left(0\right)=a,X\left(1\right)=b\end{array}$$

(7)

The recursive equations deduced from Equations (6) and (7) for different values of k is obtained as follows:

$$\begin{array}{ccccc}k=0:\hfill & & \hfill 2X\left(2\right)+\eta b+\alpha a^3+{\omega }^{2}a=0& & \hfill \end{array}$$

(8)

$$\begin{array}{ccccc}k=1:\hfill & & \hfill 6X\left(3\right)+2\eta X\left(2\right)+({\omega }^2+3\alpha a^2)b-A\Omega =0& & \hfill \end{array}$$

(9)

$$\begin{array}{cccc}k=2:\hfill & & \hfill 12X\left(4\right)+{\omega }^{2}X\left(2\right)+3\eta X\left(3\right)+3\alpha a(b^2+aX\left(2\right))=0& \end{array}$$

(10)

and so on.

The recursive relation in Equations (8)–(10) can be solved successively and then by taking the inverse differential transform $x\left(t\right)$ is obtained.

3.1. Example 1: Free Duffing Oscillator under Damping Effect

Let

$$a=1.0,b=0.0,\omega =1,\eta =0.05,\alpha =0.15,A=0$$

(11)

The analytic expansion $x\left(t\right)$ of Equation (1a) for the given values in (11) is given as follows:

$$\begin{array}{c}\hfill x\left(t\right)=1.0-0.575t^2+0.00958333t^3+0.0693594t^4\\ \hfill -0.00138839t^5-0.00830017t^6\\ \hfill +0.0002253t^7+0.00136295t^8+\cdots \end{array}$$

(12)

Figure 1 shows the comparison between the results obtained using the DTM, Equation (12), and the numerical results obtained by Runge-Kutta fourth-order accuracy method. It is clear that the results using the DTM have a reasonable agreement with the results obtained using only the fourth-order Runge-Kutta numerical method in a small range of the solution domain.

Now we are improving the accuracy of the differential transform solution using the MDTM [22]. We first apply Laplace transform to the series solution given by Equation (12) to obtain

$$\begin{array}{c}\hfill L\left[x\left(t\right)\right]=\frac{54.954}{s^9}+\frac{1.1355}{s^8}-\frac{5.9761}{s^7}-\frac{0.1666}{s^6}\\ \hfill +\frac{1.6646}{s^5}+\frac{0.0575}{s^4}-\frac{1.15}{s^3}+\frac{1}{s}\end{array}$$

(13)

As the first step of the procedure of the MDTM [22] replacing s by $1/t$, calculating the Padé approximant of [4/4] and letting $t=1/s$ gives the following:

$$\begin{array}{c}\hfill \left[\frac{4}{4}\right]=\frac{0.5525+10.6558s+0.16667s^2+s^3}{11.9025+0.6867s+11.8058s^2+0.16667s^3+s^4}\end{array}$$

(14)

Taking the inverse Laplace transform to the Padé approximant [4/4], Equation (14), to obtain the solution by the MDTM as follows:

$$\begin{array}{c}\hfill x\left(t\right)=0.9962e^{-0.0262t}cos(1.0551t)+0.0038e^{-0.0572t}cos(3.2685t)\end{array}$$

(15)

Figure 2 depicts the comparison between the MDTM results obtained by the Padé approximant of [4/4] and the results obtained using the fourth-order Runge-Kutta numerical method. It is clear that the MDTM result obtained by the real part of Padé approximate gives an excellent agreement with the result obtained using the fourth-order Runge-Kutta numerical method.

3.2. Example 2: Forced Nonlinear Duffing Oscillator with Damping Effect

In Equation (1a) let

$$\begin{array}{c}\hfill a=1.0,b=0.0,\eta =0.03,\alpha =0.03,A=0.15,\omega =1.0,\text{and}\Omega =0.8\end{array}$$

(16)

In this case the analytic expansion takes the form:

$$\begin{array}{c}\hfill x\left(t\right)=1.0-0.515t^2+0.02515t^3+0.0440155t^4-0.00219932t^5\\ \hfill -0.0015002t^6+0.0000701175t^7+0.0000273301t^8+\cdots \end{array}$$

(17)

Figure 3 illustrates the comparison of results obtained by the differential transform method which given by Equation (17) and the fourth-order Runge-Kutta numerical method. Clearly, the weakness of accuracy of the DTM result even for short time domain.

To apply the modified differential transform method, we first get the Laplace transform to the time series solution which is given by Equation (17), yields

$$\begin{array}{c}\hfill L\left[x\left(t\right)\right]=\frac{1}{s}-\frac{1.03}{s^3}+\frac{0.1509}{s^4}+\frac{1.0564}{s^5}-\frac{0.2639}{s^6}\\ \hfill -\frac{1.08015}{s^7}+\frac{0.3534}{s^8}-\frac{1.10195}{s^9}\end{array}$$

(18)

Using the concept of the MDTM [22]. Replace s by $1/t$ in Equation (18), calculate the Padé approximant of [4/4] and after that let $t=1/s$ which gives the following:

$$\begin{array}{c}\hfill \left[\frac{4}{4}\right]=\frac{0.1392+0.64s+0.03s^2+s^3}{0.6592+0.0192s+1.67s^2+0.03s^3+s^4}\end{array}$$

(19)

Applying the inverse Laplace transform to the Padé approximant of [4/4] in Equation (19) to obtain the solution by the MDTM in the form:

$$\begin{array}{c}\hfill x\left(t\right)=-0.0236cos(0.8t)+0.3832sin(0.8t)\\ \hfill +e^{-0.015t}\left(1.0236cos(1.0148t)-0.287sin(1.0148t)\right)\end{array}$$

(20)

Figure 4 shows the comparison between the MDTM result obtained by the Padé approximant of [4/4] and the result obtained using the fourth-order Runge-Kutta numerical method. The MDTM results obtained by the real part of the Padé approximant of [4/4] are clearly in excellent agreement with the result obtained using Runge-Kutta fourth-order accuracy numerical method.

4. Forced van der Pol Oscillator under Damping Effect

Nayfeh and Mook from [12] derived the Rayleigh equation in the form

$$y^{\prime \prime }+ϵ[{\displaystyle \frac{1}{3}}{\left(y^{\prime }\right)}^2-1]y^{\prime }+y=0,\text{with}{y}^{\prime }=dy/d\tau$$

(21)

Setting $x=y^{\prime }$ we find the van der Pol equation

$$x^{\prime \prime }+ϵ\left[x^2-1\right]x+x=0\text{with}x=dx/d\tau$$

(22)

which arises in the study of circuits containing vacuum tubes and given by Equation (1b) with $A=0$.

Electrical circuit involving a triode results a forced van der Pol oscillator Equation (1b) where $A\ne 0$, see Figure 5. The circuit contains: a triode, a resistor R, a coupled inductor L and mutual inductance M. In the serial $RCL$ circuit there is a current i, and towards the triode anode (plate) a current $i_a$, while there is a voltage $u_g$ on the triode control grid, see Figure 5. The van der Pol oscillator is forced by an $AC$ voltage source $E_s$.

To investigate a semi-analytic solution of forced van der Pol oscillator under damping effect we apply the MDTM to Equation (1b) with initial conditions $x\left(0\right)=a$ and $\dot{x}=b$. The differential transform of this equation has the recurrence relation

$$\begin{array}{c}\hfill (k+2)(k+1)X(k+2)+{\omega }^{2}X\left(k\right)+ϵ[(k+1)X(k+1)\\ \hfill -\sum _{n=0}^k\sum _{m=0}^n(k-n+1)X\left(m\right)X(n-m)X(k-n+1)]\\ \hfill -\frac{A{\Omega }^k}{k!}sin\left(\frac{k\pi }{2}\right)=0\end{array}$$

(23)

$$\begin{array}{c}\hfill X\left(0\right)=a,X\left(1\right)=b\end{array}$$

(24)

The recursive equations deduced from Equations (23) and (24) for $k=0,1,2$ are obtained as follows:

$$\begin{array}{ccccc}\hfill k=0:& & \hfill 2X\left(2\right)+{\omega }^{2}a+ϵ(a^2-1)b=0& & \end{array}$$

(25)

$$\begin{array}{ccccc}\hfill k=1:& & \hfill 6X\left(3\right)+{\omega }^{2}b+2ϵX\left(2\right)(a^2-1)+2ϵa{b}^2-A\Omega =0& & \end{array}$$

(26)

$$\begin{array}{cccc}\hfill k=2:& & \hfill 12X\left(4\right)+{\omega }^{2}X\left(2\right)+ϵ[3(a^2-1)X\left(3\right)+6abX\left(2\right)+b^3]=0& \end{array}$$

(27)

and so on.

4.1. Example 3:

Let

$$\begin{array}{c}\hfill a=1.0,b=0.0,\omega =1.0,ϵ=0.04,A=0.04,\text{and}\Omega =1.4\end{array}$$

(28)

For the values given in Equation (28) the analytic expansion of the nonlinear van der Pol Equation (1b) is given by:

$$\begin{array}{c}\hfill x\left(t\right)=1.0-0.5t^2+0.009333t^3+0.041667t^4{t}^4-0.003381t^5\\ \hfill -0.001327t^6+0.000599t^7-0.000009t^8+\cdots \end{array}$$

(29)

To use the modified differential transform method we take the Laplace transform of the series solution in Equation (29), yields:

$$\begin{array}{c}\hfill L\left[x\left(t\right)\right]=\frac{1}{s}-\frac{1}{s^3}+\frac{0.056}{s^4}+\frac{1}{s^5}-\frac{0.40576}{s^6}-\frac{0.9552}{s^7}+\frac{3.01838}{s^8}\\ \hfill -\frac{0.353677}{s^9}-\frac{25.0076}{s^{10}}+\frac{31.4072}{s^{11}}\end{array}$$

(30)

Forthwith as in [22] we replace s by $1/t$, calculating the Padé approximant of [3/3] and [4/4] and letting $t=1/s$ gives us the following:

$$\begin{array}{c}\hfill \left[\frac{3}{3}\right]=\frac{5.24571+93.6735s+s^2}{93.6175+6.24571s+93.6735s^2+s^3}\end{array}$$

(31)

$$\begin{array}{c}\hfill \left[\frac{4}{4}\right]=\frac{0.13895+7.72697s+0.168318s^2+s^3}{7.71754+0.251268s+8.72697s^2+0.168318s^3+s^4}\end{array}$$

(32)

Taking the inverse Laplace transform of the Padé approximant of [3/3] and [4/4] in Equations (31) and (32) we obtain the semi-analytic solutions respectively, as follows:

$$\begin{array}{c}\hfill x\left(t\right)\approx e^{-0.02802t}(cost+0.02803sint),\\ \hfill x\left(t\right)\approx e^{-0.0062t}(cost+0.01264sint)\end{array}$$

(33)

$$\begin{array}{c}\hfill +e^{-0.07798t}(0.0003cos(2.7785t)-0.0023sin(2.7785t))\end{array}$$

(34)

Figure 6a illustrated to show a comparison of the solution by the differential transform method, Equation (29), and the numerical solution by the fourth-order Runge-Kutta method. It is clear that, the result obtained by DTM have not reasonable agreement with the numerical result by Runge-Kutta for a long time domain.

Figure 6b shows the comparison between the MDTM results obtained by the real part of Padé approximant of orders [3/3] and [4/4] whose given by Equations (33) and (34), respectively, and the result obtained by the fourth-order Runge-Kutta numerical method. The MDTM result obtained by the real part of Padé approximant [3/3] in Equation (33) shows some discrepancies in comparison to the result obtained using the fourth-order Runge-Kuta numerical method. However, it is clear that the result of the result of the MDTM by the real part of Padé approximant [4/4] in Equation (34) has excellent agreement and seems to coincide with the numerical result.

4.2. Example 4:

Let

$$\begin{array}{c}\hfill a=1.0,b=0.0,ϵ=0.004,A=0.9,\omega =2.0\text{and}\Omega =1.4\end{array}$$

(35)

The nonlinear differential Equation (1b) in this case has the analytic solution in the form:

$$\begin{array}{c}\hfill x\left(t\right)=1-2t^2+0.21t^3+0.666667t^4-0.06578t^5-0.088329`t^6\\ \hfill +0.010248t^7+0.005618t^8-0.002084t^9+0.000145t^{10}\\ \hfill +0.000545t^{11}-0.000143t^{12}+\cdots \end{array}$$

(36)

To apply the MDTM we first take Laplace transform of the DTM solution given in Equation (36), yields:

$$\begin{array}{cc}\hfill L\left[x\left(t\right)\right]=& \frac{1}{s}-\frac{4}{s^2}+\frac{1.26}{s^4}+\frac{16}{s^5}-\frac{7.8936}{s^6}-\frac{63.5968}{s^7}+\frac{51.647808}{s^8}\hfill \\ & +\frac{226.503501}{s^9}-\frac{756.140505}{s^{10}}+\frac{524.935235}{s^{11}}\hfill \\ & +\frac{21755.582199}{s^{12}}-\frac{68682.662616}{s^{13}}\hfill \end{array}$$

(37)

Replace s by $1/t$ in Equation (37) and calculate the Padé approximants [4/4] and [6/6] and after that let $t=1/s$ to obtain

$$\begin{array}{c}\hfill \left[\frac{4}{4}\right]=\frac{0.554007+3.493991s-0.860981s^2-s}{15.0608-1.629917s-0.506009s^2-0.860981s^3-s^4}\end{array}$$

(38)

$$\begin{array}{c}\hfill \left[\frac{6}{6}\right]=\frac{45.6864+71.9781s+3.2114s^2+38.5489s^3+0.7508s^4+s^5}{283.4087+7.7739s+225.2278s^2+4.9547s^3+42.5489s^4+0.7508s^5+s^6}\end{array}$$

(39)

The inverse Laplace transform of (38) and (39), respectively, gives the semi-analytical solutions as follows:

$$\begin{array}{c}\hfill x\left(t\right)=-0.0058e^{-2.2603t}+0.01583e^{1.6398t}\\ \hfill +e^{-0.12024t}(0.9899cos(2.0122t)+0.0394sin(2.0122t))\end{array}$$

(40)

$$\begin{array}{c}\hfill x\left(t\right)=e^{-0.00698t}(0.0091cos(1.3895t)+0.4314sin(1.3895t))\\ \hfill +e^{-0.00015t}(cos(2.00267t)-0.29905sin92.00267t))\end{array}$$

(41)

Figure 7 illustrates the comparison of the fourth-order Runge-Kutta numerical solution and the solution of the MDTM by the real part of the Padé approximant of order [4/4]. In Figure 8, it is clear that the solution by the MDTM is unstable. For this reason we have another attempt to obtain a more accurate and stable semi-analytical solution. To resolve this problem we use again the MDTM but with the real part of the Padé approximant of order [6/6] given by Equation (41).

Example 4 shows that the solution by the MDTM with the real part of the Padé approximant [6/6] not only matches perfectly with the numerical solution for a long time domain, as in Figure 8, but also shows that the phase plane given by the two methods seems to be identical, show Figure 9.

5. Conclusions

The main goal of researchers who are interested in solving nonlinear differential equations is to obtain analytical solutions along with numerical solutions. These researchers have relied on some methods such as the multiple time scales method and the harmonic balance method and others. Other researchers have taken another turn and used the modified differential to transform method (MDTM) to obtain semi-analytic solutions of free non-linear oscillation by adding Laplace transform and Padé approximant [4/4]. Here, we extend their studies and provided semi-analytical solutions of forced oscillations of Duffing and van der Pol under damping effects.

Through some applications, we have provided records to indicate the success of the semi-analytical solution by the MDTM of nonlinear differential equations. The study focused on van der Pol and Duffing nonlinear oscillators under damping effects, due to their importance in science and engineering.

The numerical results demonstrate the validity and applicability of analytic solutions that we obtained by using the MDTM with appropriate values of parameters. This assures us the extent of success of using the modified differential to transform method to obtain analytical solutions of non-linear differential equations.

To get the analytical solution using the MDTM method one can obtain the result of the real part of Padé approximant of any order like [3/3], [4/4], [5/5] and higher where at least one of these approximant gives a satisfactory accurate analytic solution.