On the Multistage Differential Transformation Method for Analyzing Damping Duffing Oscillator and Its Applications to Plasma Physics

Abstract

The multistage differential transformation method (MSDTM) is used to find an approximate solution to the forced damping Duffing equation (FDDE). In this paper, we prove that the MSDTM can predict the solution in the long domain as compared to differential transformation method (DTM) and more accurately than the modified differential transformation method (MDTM). In addition, the maximum residual errors for DTM and its modification methods (MSDTM and MDTM) are estimated. As a real application to the obtained solution, we investigate the oscillations in a complex unmagnetized plasma. To do that, the fluid govern equations of plasma species is reduced to the modified Korteweg–de Vries–Burgers (mKdVB) equation. After that, by using a suitable transformation, the mKdVB equation is transformed into the forced damping Duffing equation.

1. Introduction

Mathematical techniques are very important tools in mathematics. Mathematicians have developed many mathematical methods to compute linear or nonlinear differential equations which describe many important phenomena and applications in science [1,2,3,4,5,6,7]. The mathematical techniques are classified as algebraic methods, semi-approximate, general analytical, approximate analytical, numerical, or qualitative techniques. The basic concept of approximate analytical techniques such as Adomian decomposition method (ADM), Laplacian decomposition method (LDM), or differential transformation method (DTM) is assuming that the solution is descried by a Taylor expansion form. Indeed, some solutions of equations have well-known Taylor expansions such as exponential function or hyperbolic function. In this case, it is easy to determine the exact solution by a few terms of the Taylor expansion series. Otherwise, the approximate solution will be obtained in the form of few terms of Taylor expansion series. Since Taylor expansion is local convergent about the initial condition, the method can approximate the solution in the neighborhood of the initial point. Thus, the solution is obtained in a very short domain. This feature of ADM, LDM or DTM has been mentioned by some researchers [8,9,10,11,12]. DTM has been improved by dividing the domain into subdomains and modifying the initial point in each subdomain. The other modification is by using the Laplacian transformation and Padé approximate. In Section 2, we describe these modifications in details. However, it is very important to determine the optimal modification of DTM to present fast and accurate techniques.

Some of the most important and famous differential equations whose solutions are related to many natural phenomena, physical concepts, and engineering phenomena are the Duffing equation (including conservative and non-conservative cases), the Helmholtz equation (including conservative and non-conservative cases), and their families [13,14,15,16,17,18,19,20,21,22,23,24]. Given the importance of the family of Duffing equation, a great effort has been made by many researchers to solve this equation and its family, with a numerical, analytical, or semi-analytical solution according to the type of Duffing equation. Examples of these approximate methods for solving the conservative Duffing equation ($u^{″}+\beta u\left(t\right)+\gamma u^3=0$) include the homotopy perturbation method [25], harmonic balance method [26], energy balance method [27], modified variational approach [28], and coupled homotopy–variational approach [29]. On the other hand, many researchers have tried to find a solution to the damping Duffing equation (DDE), $\left(u^{\prime }+\alpha u^{\prime }+\beta u\left(t\right)+\gamma u^3=0\right)$ [30,31,32,33,34,35], since it is more closely related to reality than the undamping Duffing equation, which is correct only for idealized isolated systems, i.e., systems in which the frictional force and viscosity are absent. One of the most important approximate methods that has been used and developed to solve many differential equations is DTM, which has been used in solving DDE [11]. Nourazar et al. [11] used the modified DTM to get an approximate solution to the DDE. The authors compared their solution with both the fourth-order Runge–Kutta (RK4) numerical solution and the DTM solution. They found that the DTM solution is suitable only for small time intervals while the MDTM solution is suitable for the whole time domain. In our study, we solve the forced damping Duffing equation (FDDE) $\left(u^{″}+\alpha u^{\prime }+\beta u\left(t\right)+\gamma u^3=F\right)$ using the multistage differential transformation method (MSDTM) for arbitrary initial conditions. Moreover, we compare the approximate solutions of DTM and MDTM as well as the numerical solution using RK4 in order to determine the optimal technique. Furthermore, the oscillations in complex unmagnetized plasmas are investigated by reducing the fluid govern equations of the plasma species to an evolution equation and then transform this equation to the Duffing-type equation using a suitable transformation.

2. Methodology

This section is devoted to briefly describing DTM and its modifications. Assume the following ordinary differential equation (ODE)

$$P(u,u^{\prime },u^{″},\dots .)=0,$$

(1)

where $u\left(t\right)$ is the solution of this ODE in domain $[t_0,t_N]$, P is a polynomial in terms of u and its derivative, and $u\left(t_0\right)=c$.

2.1. Differential Transformation Method (DTM)

Assume that the goal is finding the approximate solution of Equation (1). The main concept of DTM is based on applying the differential transformation $u\left(t\right)⟹U\left(k\right)$ at $t=t_0$ as follows:

$$U\left(k\right)=\frac{1}{k!}{\left[\frac{d^{k}u\left(t\right)}{d{t}^k}\right]}_{t=t_0}.$$

(2)

The differential inverse transformation $U\left(k\right)⟹u\left(t\right)$ is defined as

$$u\left(t\right)=\sum _{k=0}^{\infty }U\left(k\right){(t-t_0)}^k,$$

(3)

Inserting Equation (2) into Equation (3), $u\left(t\right)$ can be approximated in finite number series as follows

$$u_N\left(t\right)=\sum _{k=0}^N\frac{{(t-t_0)}^k}{k!}{\left[\frac{d^{k}u\left(t\right)}{d{t}^k}\right]}_{t=t_0}=g_N.$$

(4)

Some differential transformation rules are introduced in

Table 1. Differential transformation rules.

Original Function	Transformed Function
$u\left(t\right)\pm v\left(t\right)$	$U\left(t\right)\pm V\left(t\right)$
$cu\left(t\right)$	$cU\left(t\right)$ (c is constant)
$\frac{du\left(t\right)}{dt}$	$(k+1)U(k+1)$
$\frac{d^{n}u\left(t\right)}{d{t}^n}$	$\frac{(k+n)!}{k!}U(k+n)$
$u\left(t\right)v\left(t\right)$	${\sum }_{m=0}^{k}U\left(m\right)V(k-m)$

.

It is well known that, since the DTM based on Taylor expansion, the approximate solution if it is locally analytic converges to the exact solution with the following approximated error

$${|u\left(t\right)}_{-}{g}_N\left(t\right)|\le \frac{M}{(N+1)!}{|t-t_0|}^{N+1},$$

where $\left|u\right(t_N\left)\right|\le M$.

It is obvious that the error increases when $|t-t_0|$ incenses for fixed term N. Note that DTM gives accurate results only in a small domain around the initial point. Therefore, to obtain good results, some modifications to this method must be introduced. There are some attempts to improve this method, such as the modified differential transformation method (MDTM) and the multistage differential transformation method (MSDTM).

2.2. Modified Differential Transformation Method (MDTM)

MDTM is presented in [11]. The idea is described simply by applying the Laplacian transformation into Equation (3) $\mathcal{L}u\left(t\right)$. We obtain the polynomial in terms of $1/t^s$. Next, we use Padé approximate, $[3/3]$ or $[4/4]$, and then apply the Laplacian inverse transform. The method is improved and able to approximate the solution in long domain.

Definition 1.

We say the function $g\left(t\right)$ is Padé approximate of order $[m/n]$ for function $u\left(t\right)$ if

$$g\left(t\right)=\frac{a_0+a_{1}t+a_2{t}^2+\dots \dots +a_m{t}^m}{1+b_0+b_{1}t+b_2{t}^2+\dots \dots +b_n{t}^n},$$

where $u\left(0\right)=g\left(0\right),u^{\prime }\left(0\right)=g^{\prime }\left(0\right),u^{″}\left(0\right)=g^{″}\left(0\right),\dots \dots ,u^{(m+n)}\left(0\right)=g^{(m+n)}\left(0\right)$. The constants $a_i,i=1,2,\dots ,m$ and $b_j,j=1,2,\dots ,n$ are uniquely determined. The Padé approximate is unique for given n and m.

2.3. Multistage Differential Transformation Method (MSDTM)

The other modification is MSDTM. The main concept is dividing the domain into subdomains $[t_i,t_{i+1}]=D_i$ and applying DTM in each subdomain with the initial condition at $t_i$ to approximate $u\left(t\right)$ at the subdomain $D_i$.

2.4. Example

In this section, we apply DTM and its modifications to one of the most famous equations in dynamic systems which is called the Duffing oscillator or the Duffing equation. It is known that the Duffing equation has many formulas, and, in this paper, we restrict our attention to investigating the forced damping Duffing equation (FDDE) $\left(u^{″}+\alpha u^{\prime }+\beta u\left(t\right)+\gamma u^3=F\right).$ This equation is non-integrable and does not have an exact solution except under certain conditions on its coefficients $\left(\alpha ,\beta ,\gamma \right)$. Therefore, the approximate solution to the following FDDE for arbitrary values of its coefficients $\left(\alpha ,\beta ,\gamma \right)$ and for arbitrary initial conditions is obtained:

$$\left\{\begin{array}{c}{u}^{″}+\alpha u^{\prime }+\beta u\left(t\right)+\gamma u^3=F,\\ u\left(0\right)=u_0\&u^{\prime }=u_0^{\prime }.\end{array}\right.$$

(5)

In the following analysis, we give some numerical examples to solve the initial value problem (i.v.p.) (5) using the aforementioned methods and examine the accuracy of these methods for calculating the residual error for each method compared to the RK4.

2.4.1. MDTM

Firstly, let us use the same values of $\left(\alpha ,\beta ,\gamma ,F\right)=\left(0.5,25,25,0\right)$ as mentioned by [11] with the initial conditions $u\left(0\right)=0.1$ and $u^{\prime }\left(0\right)=0.$ Note that the solution of the i.v.p. (5) for unforced ($F=0$) using MDTM is introduced in details in [11]. In the case of using Padé approximate of [3/3], we have

$$\begin{array}{cc}\hfill u\left(t\right)& =0.00194+0.000238e^{-0.25t}(411cos\left(5.068t\right)\hfill \\ \hfill & +20.273sin\left(5.068t\right)),\hfill \end{array}$$

(6)

In the case of using Padé approximate of [4/4], the solution is approximated as

$$\begin{array}{cc}\hfill u\left(t\right)& =A{e}^{(-0.60107-15.0816i)t}+B{e}^{(-0.60107+15.0816i)t}\hfill \\ \hfill & +C{e}^{(-0.24894-5.0125i)t}+D{e}^{(-0.24894+5.0125i)t},\hfill \end{array}$$

(7)

with

$$\begin{array}{cc}\hfill A& =1.6932\times {10}^{-5}-1.3567\times {10}^{-4}i,\hfill \\ \hfill B& =1.6932\times {10}^{-5}+1.3567\times {10}^{-4}i,\hfill \\ \hfill C& =4.9983\times {10}^{-2}-2.5633\times {10}^{-3}i,\hfill \\ \hfill D& =4.9983\times {10}^{-2}+2.5633\times {10}^{-3}i.\hfill \end{array}$$

In the second example, we use the values $\left(\alpha ,\beta ,\gamma ,F\right)=\left(1,20,2,0\right)$ and with initial condition $u\left(0\right)=-0.2$ and $u^{\prime }\left(0\right)=2$ and apply Padé approximate of [3/3] and $[4/4]$. The solution in the case of using Padé approximate of [3/3] reads

$$\begin{array}{cc}\hfill u\left(t\right)& =0.003101exp(-6.3493t)\hfill \\ \hfill & +exp(-0.52098t)(0.434516sin(4.4046t)\hfill \\ \hfill & -0.203101cos\left(4.4046t\right)),\hfill \end{array}$$

(8)

and for $[4/4]$ reads

$$\begin{array}{cc}\hfill u\left(t\right)& =A{e}^{(-2.0169+12.6572i)t}+B{e}^{(-2.0169-12.6572i)t}\hfill \\ \hfill & +C{e}^{(-0.4965+4.4826i)t}+D{e}^{(-0.4965-4.4826i)t},\hfill \end{array}$$

(9)

with

$$\begin{array}{cc}\hfill A& =2.265\times {10}^{-4},-4.807\times {10}^{-5}i,\hfill \\ \hfill B& =2.265\times {10}^{-4},+4.807\times {10}^{-5}i,\hfill \\ \hfill C& =-0.100226-0.21195i,\hfill \\ \hfill D& =-0.100226+0.21195i.\hfill \end{array}$$

2.4.2. MSDTM

In this work, we focus our attention to solve the i.v.p. (5) for arbitrary initial conditions using MSDTM by dividing the domain $[0,20]$ to subdomains with time step ${10}^{-2}$ and apply DTM with $k=3$ to find $u^i$ as follows:

$$u_{k+1}^i=\frac{k!}{(k+1)!}{y}_k^i,$$

(10)

$$y_{k+1}^i=\frac{k!}{(k+1)!}\left[-\beta u_k^i-\alpha y_k^i-\gamma \sum _{r=0}^k\left(\sum _{l=0}^r\left(u_l^i{u}_{r-l}^i\right)\right)u_{k-r}^i+F\right],$$

(11)

where $y=u^{\prime }$.

To check the accuracy of the aforementioned methods as compared the RK4 solution, we use the following error formula for the maximum residual error

$$L_D\left(\mathrm{method}\right)=\underset{t_0\le t\le t_N}{max}\left|\mathrm{RK}\left(t\right)-u\left(t\right)\right|.$$

Figure 1 and Figure 2 demonstrate the approximate solutions to the i.v.p. (5) for different values of the coefficients $\left(\alpha ,\beta ,\gamma ,F\right)$. The results show that the MDTM4 and MSDTM are better approximations than MDTM3. Moreover, the comparison of the maximum residual errors for the approximate solutions shown in

Table 2. The error $L_D\left(\mathrm{methods}\right)$ is estimated for different values of the coefficients $\left(\alpha ,\beta ,\gamma ,u_0,u_0^{\prime }\right).$

$\left(\mathit{\alpha },\mathit{\beta },\mathit{\gamma },{\mathit{u}}_0,{\mathit{u}}_0^{\prime }\right)$	Time Range	${\mathit{L}}_{\mathit{D}}$(MDTM3)	${\mathit{L}}_{\mathit{D}}$(MDTM4)	${\mathit{L}}_{\mathit{D}}$(MSDTM)
$\left(0.5,25,25,0.1,0\right)$	$0\le t\le 20$	$1.19631\times {10}^{-2}$	$1.67636\times {10}^{-3}$	$4.81974\times {10}^{-4}$
$\left(1,20,2,-0.2,2\right)$	$0\le t\le 6$	$1.83182\times {10}^{-2}$	$7.99595\times {10}^{-3}$	$1.96895\times {10}^{-5}$

proves that the accurate method is MSDTM. Aljahdaly [10] proved that the MSDTM and RK4 techniques have the same accuracy, but MSDTM is faster than RK4. Thus, we conclude that MSDTM is a fast, accurate, and reliable method for many differential equations in physics and in different branches of science. In the next section, a new application to the damping Duffing equation in plasma physics is introduced.

3. Application in Plasma Physics

Let us consider the propagation of nonlinear structures in a complex unmagnetized plasma composed of inertial positive ions (with subscript “i”) and two different types of electrons (with subscripts “l” and “h” for the lower and higher electron temperatures, respectively) that follow the kappa distribution in addition to static dust grains with negative charge (with subscript “d”) [36]. Accordingly, the neutrality condition reads $n_l^{\left(0\right)}+n_h^{\left(0\right)}+z_d{n}_d^{\left(0\right)}=n_i^{\left(0\right)}$, where $n_j^{\left(0\right)}$ represents the unperturbed number density of species j $\left(j\equiv l,h,d,i\right)$ and $z_d$ gives the number of electrons residing on the surface of the dust grains. The dynamics of the nonlinear structures whose phase speed is much larger than the ion thermal speed but smaller than the electron thermal speed are governed by the following dimensionless fluid continuity, momentum, and Poisson’s equations, respectively,

$$\left\{\begin{array}{c}{\partial }_t{n}_i+{\partial }_x\left(n_i{u}_i\right)=0,\\ {\partial }_t{u}_i+u_i{\partial }_x{u}_i+{\partial }_x\varphi =\eta {\partial }_x^2{u}_i,\\ {\partial }_x^2\varphi -n_e+n_i-{\mu }_d=0,\end{array}\right.$$

(12)

where the number density of the electrons in kappa distribution is given by

$$\begin{array}{cc}\hfill n_e& =n_l+n_h={\mu }_l{\left(1-\frac{{\sigma }_l\varphi }{R_l}\right)}^{S_l}+{\mu }_h{\left(1-\frac{{\sigma }_h\varphi }{R_h}\right)}^{S_h}\hfill \\ \hfill & \equiv {\Gamma }_0+{\Gamma }_1\varphi +{\Gamma }_2{\varphi }^2+{\Gamma }_3{\varphi }^3+\cdots ,\hfill \end{array}$$

(13)

with

$$\begin{array}{cc}\hfill {\Gamma }_0& ={\mu }_l+{\mu }_h,\hfill \\ \hfill {\Gamma }_1& =-\left[\frac{S_l{\mu }_l{\sigma }_l}{R_l}+\frac{S_h{\mu }_h{\sigma }_h}{R_h}\right],\hfill \\ \hfill {\Gamma }_2& =\left[\frac{S_l{\mu }_l{\sigma }_l^2\left(S_l-1\right)}{2R_l^2}+\frac{S_h{\mu }_h{\sigma }_h^2\left(S_h-1\right)}{2R_h^2}\right],\hfill \\ \hfill {\Gamma }_3& =-\left[\frac{S_l{\mu }_l{\sigma }_l^3\left(S_l-1\right)\left(S_l-2\right)}{6R_l^3}+\frac{S_h{\mu }_h{\sigma }_h^3\left(S_h-1\right)\left(S_h-2\right)}{6R_h^3}\right],\hfill \\ \hfill S_l& =\left(-{\kappa }_l+\frac{1}{2}\right)\&S_h=\left(-{\kappa }_h+\frac{1}{2}\right),\hfill \\ \hfill R_l& =\left({\kappa }_l-\frac{3}{2}\right)\&R_h=\left({\kappa }_h-\frac{3}{2}\right).\hfill \end{array}$$

where $n_i/n_l/n_h$ is the normalized number density of the positive ions/low temperature electrons/high temperature electrons, $u_i$ refers to the normalized velocity of the positive ions, $\varphi$ is the normalized electrostatic potential, $\eta$ represents the normalized coefficient of ionic kinematic viscosity, ${\sigma }_{l,h}=T_{eff}/T_{l,h}$ is the electron temperature ratio, the effective electron temperature is $T_{eff}=n_e^{\left(0\right)}{T}_l{T}_h/\left(n_l^{\left(0\right)}{T}_h+n_h^{\left(0\right)}{T}_l\right)$, $n_e^{\left(0\right)}\equiv \left(n_l^{\left(0\right)}+n_h^{\left(0\right)}\right)$ is the total unperturbed electrons density, ${\mu }_d=z_d{n}_d^{\left(0\right)}/n_i^{\left(0\right)}$ is the dust concentration, ${\mu }_l=n_l^{\left(0\right)}/n_i^{\left(0\right)}$ is the concentration of low electron temperature, ${\mu }_h=n_h^{\left(0\right)}/n_i^{\left(0\right)}$ is the concentration of high electron temperature, and ${\kappa }_{l,h}\left(>3/2\right)$ is the kappa index parameter [36].

To model and analyze the nonlinear structures that can propagate in the present plasma model, the reductive perturbation method (RPM) [37,38] is used to reduce the basic set of fluid Equations (12) and (13) to an evolution equation. According to this method, the independent variables $\left(x,t,\eta \right)$ can be stretched as follows:

$$X=\epsilon \left(x-v_{ph}t\right),T={\epsilon }^{3}t\&\eta =\epsilon \tilde{\eta },$$

(14)

where $\epsilon$ is a real and small parameter ($\epsilon <<1$) that measures the strength of the nonlinearity and $v_{ph}$ represents the normalized phase velocity, which is scaled by $C_i$. In addition, the dependent quantities $\Pi \left(x,t\right)\equiv \left(n_i,u_i,\varphi \right)$ are expanded as follows:

$$\Pi \left(x,t\right)={\Pi }^{\left(0\right)}+\sum _{s=1}^{\infty }{\epsilon }^s{\Pi }^{\left(s\right)}\left(X,T\right),$$

(15)

where ${\Pi }^{\left(0\right)}\equiv {\left[1,0,0\right]}^T$, ${\Pi }^{\left(s\right)}\left(X,T\right)\equiv {\left[n_i^{\left(s\right)},u_i^{\left(s\right)},{\varphi }^{\left(s\right)}\right]}^T$, and T gives the matrix transpose.

Inserting both stretching (14) and expansion (15) into the basic set of fluid Equations (12) and (13), we get a system of reduced equations in different powers of $\epsilon$. From the lowest-order of $\epsilon$, i.e., $O\left(\epsilon \right)$, the values of the first-order quantities $\left(n_i^{\left(1\right)},u_i^{\left(1\right)}\right)$ and the phase velocity $v_{ph}$ can be obtained as

$$\begin{array}{cc}\hfill u_i^{\left(1\right)}& =v_{ph}{n}_i^{\left(1\right)}=\frac{1}{v_{ph}}{\varphi }^{\left(1\right)},\hfill \\ \hfill v_{ph}& =\frac{1}{\sqrt{{\Gamma }_1}}.\hfill \end{array}$$

(16)

The solution of next-order of $\epsilon ,$ i.e., $O\left({\epsilon }^2\right)$, gives the values of the second-order quantities $\left(n_i^{\left(2\right)},u_i^{\left(2\right)}\right)$

$$\begin{array}{cc}\hfill n_i^{\left(2\right)}& =\frac{1}{v_{ph}^4}\left(\frac{3}{2}{\varphi }^{\left(1\right)2}+v_{ph}^2{\varphi }^{\left(2\right)}\right),\hfill \\ \hfill u_i^{\left(2\right)}& =\frac{1}{v_{ph}^3}\left(\frac{1}{2}{\varphi }^{\left(1\right)2}+v_{ph}^2{\varphi }^{\left(2\right)}\right),\hfill \end{array}$$

(17)

and the Poisson’s equation gives

$$A{\varphi }^{\left(1\right)2}+A_c{\varphi }^{\left(2\right)}=0,$$

(18)

where $A=\left[3/\left(2v_{ph}^4\right)-{\Gamma }_2\right]=0$ at the critical value of low electron temperature concentration ${\mu }_l=$${\mu }_{lc}$ and the coefficient $A_c=\left(1/v_{ph}^2-{\Gamma }_1\right)$ represents the compatibility condition, i.e., $A_c=0$.

From the next-order of $\epsilon$, i.e., $O\left({\epsilon }^3\right)$, we get

$$\begin{array}{cc}\hfill {\partial }_T{n}_i^{\left(1\right)}+{\partial }_X\left(n_i^{\left(1\right)}{u}_i^{\left(2\right)}\right)+{\partial }_X\left(n_i^{\left(2\right)}{u}_i^{\left(1\right)}\right)-v_{ph}{\partial }_X{n}_i^{\left(3\right)}+{\partial }_X{u}_i^{\left(3\right)}& =0,\hfill \end{array}$$

(19)

$$\begin{array}{cc}\hfill {\partial }_T{u}_i^{\left(1\right)}+{\partial }_X\left(u_i^{\left(1\right)}{u}_i^{\left(2\right)}\right)+{\partial }_X\left(n_i^{\left(2\right)}{u}_i^{\left(1\right)}\right)-v_{ph}{\partial }_X{u}_i^{\left(3\right)}+{\partial }_X{\varphi }^{\left(3\right)}-\tilde{\eta }{\partial }_X^2{u}_i^{\left(1\right)}& =0,\hfill \end{array}$$

(20)

and the Poisson’s equation gives

$${\partial }_X\left(n_i^{\left(3\right)}-{\Gamma }_3{\varphi }^{\left(1\right)3}-2{\Gamma }_2{\varphi }^{\left(1\right)}{\varphi }^{\left(2\right)}-{\Gamma }_1{\varphi }^{\left(3\right)}+{\partial }_X^2{\varphi }^{\left(1\right)}\right)=0.$$

(21)

By solving Equations (19)–(21) with the help of Equations (16) and (17), we finally get the mKdVB equation

$${\partial }_T\phi +P_1{\phi }^2{\partial }_X\phi +P_2{\partial }_X^3\phi =P_3{\partial }_X^2\phi ,$$

(22)

with

$$\begin{array}{cc}\hfill P_1& =\left(15-6{\Gamma }_3{v}_{ph}^6\right)/\left(4v_{ph}^3\right),\hfill \\ \hfill P_2& =\frac{v_{ph}^3}{2}\&P_3=\frac{\tilde{\eta }}{2},\hfill \end{array}$$

where $\phi \equiv {\varphi }^{\left(1\right)}.$

It is known that Equation (22) supports the shock solution due to the presence of ion kinematic viscosity. However, in this paper, we want to investigate the damping oscillations in the present model. Accordingly, the transformation $\phi \left(X,T\right)=\mathsf{\Phi }\left(\xi \right)$, where $\xi =X+v_{f}T,$ is used to transform Equation (22) into the FDDE as follows:

$${\phi }^{″}+\alpha {\phi }^{\prime }+\beta \phi +\gamma {\phi }^3=F,$$

(23)

where $\alpha =-P_3/P_2$, $\beta =v_f/P_2$, $\gamma =P_1/\left(3P_2\right)$, and F is the constant of integration.

Let us now investigate the effect of typical complex plasma parameters, namely $\left({\kappa }_l,{\sigma }_h,{\mu }_l,{\mu }_h\right)=\left(3,0.1,{\mu }_c,0.4\right)$, and different values for $\left({\sigma }_l,{\kappa }_h,\tilde{\eta }\right)$ on the profile of plasma oscillations. Some plasma data are used as an example for investigating the solution of MSDTM, as shown in Figure 3. It is clear from the results in Figure 3 that the enhancement of the viscosity parameter $\tilde{\eta }$ leads to an increase in the number of oscillations and decreasing the time of periodicity. Note that the effect ${\sigma }_l$ has on the profile of oscillation is the same as its effect on $\tilde{\eta }$ while ${\kappa }_h$ has the opposite effect, in which the number of oscillations decreases and the time periodicity increases with the enhancement of ${\kappa }_h$.

4. Conclusions

The forced damping Duffing equation $\left({\phi }^{″}+\alpha {\phi }^{\prime }+\beta \phi +\gamma {\phi }^3=F\right)$ with arbitrary initial conditions is investigated numerically via the highly-accurate MSDTM. The comparison between the approximate solutions using MDTM and MSDTM with RK4 numerical solution is reported. Moreover, the maximum residual error for all approximate numerical solutions as compared to the RK4 solution is estimated. It is observed that the approximate numerical solution using MSDTM is highly accurate and better than both DTM and MDTM. Furthermore, the application of the FDDE in the practical plasma model is investigated to study the dynamics of nonlinear oscillations that occur in a complex unmagnetized plasma. This solution might help many researchers in studying and investigating many problems in various fields of science such as plasma physics and optical fiber.

Future work: in this work, the MSDTM is devoted for solving the FDDE for constant force, but sometimes the perturbation force is not constant but periodic with time $\left({\phi }^{″}+\alpha {\phi }^{\prime }+\beta \phi +\gamma {\phi }^3=f\left(t\right)\right)$, this is considered an important and vital problem but out of the present scope.