Two Reliable Methods for The Solution of Fractional Coupled Burgers’ Equation Arising as a Model of Polydispersive Sedimentation

Taking account into the following nonlinear time fractional equation in order to explain the basic idea of the implemented method [15]

Fu,∂<!-- ? -->α<!-- a -->u∂<!-- ? -->tα<!-- a -->,∂<!-- ? -->u∂<!-- ? -->x,∂<!-- ? -->u∂<!-- ? -->y,∂<!-- ? -->2α<!-- a -->u∂<!-- ? -->t2α<!-- a -->,∂<!-- ? -->2u∂<!-- ? -->x2,∂<!-- ? -->2u∂<!-- ? -->y2,…<!-- ? -->=0 $$\begin{array}{} \displaystyle F\left( u,\frac{\partial ^{\alpha }u}{\partial t^{\alpha }},\frac{\partial u% }{\partial x},\frac{\partial u}{\partial y},\frac{\partial ^{2\alpha }u}{% \partial t^{2\alpha }},\frac{\partial ^{2}u}{\partial x^{2}},\frac{\partial ^{2}u}{\partial y^{2}},\ldots \right) =0 \end{array}$$(6)

where the fractional derivatives are in conformable sense. We can introduce the wave variable as

u(x,y,t)=u(η<!-- ? -->),η<!-- ? -->=kx+wy+ctα<!-- a -->α<!-- a --> $$\begin{array}{} \displaystyle u(x,y,t)=u(\eta ),\eta =kx+wy+c\frac{t^{\alpha}}{\alpha} \end{array}$$(7)

where k, w, c are arbitrary constants that can be examined later. With the help of conformable chain rule [1], we have

∂<!-- ? -->α<!-- a -->(.)∂<!-- ? -->tα<!-- a -->=cd(.)dη<!-- ? -->,∂<!-- ? -->(.)∂<!-- ? -->x=kd(.)dη<!-- ? -->,∂<!-- ? -->(.)∂<!-- ? -->y=wd(.)dη<!-- ? -->,…<!-- ? -->. $$\begin{array}{} \displaystyle \frac{\partial ^{\alpha}(.)}{\partial t^{\alpha }}=c\frac{d(.)}{d\eta },~% \frac{\partial (.)}{\partial x}=k\frac{d(.)}{d\eta },~\frac{\partial (.)}{% \partial y}=w\frac{d(.)}{d\eta },\ldots . \end{array}$$(8)

Hence Eq.(6) changes into differential equation with integer order as follows.

Qu,uη<!-- ? -->,uη<!-- ? -->η<!-- ? -->,uη<!-- ? -->η<!-- ? -->η<!-- ? -->,uη<!-- ? -->η<!-- ? -->η<!-- ? -->η<!-- ? -->,...=0. $$\begin{array}{} \displaystyle Q\left( u,u_{\eta },u_{\eta \eta },u_{\eta \eta \eta },u_{\eta \eta \eta \eta },...\right) =0. \end{array}$$(9)

Due to exp-function method, it is supposed that the wave solution can be regarded in the following form

u(η<!-- ? -->)=∑<!-- ? -->n=−<!-- - -->djanenη<!-- ? -->∑<!-- ? -->s=−<!-- - -->qpbmemη<!-- ? -->=ajejη<!-- ? -->+…<!-- ? -->+a−<!-- - -->de−<!-- - -->dη<!-- ? -->bpepη<!-- ? -->+…<!-- ? -->+a−<!-- - -->qe−<!-- - -->qη<!-- ? --> $$\begin{array}{} \displaystyle u(\eta )=\frac{\sum_{n=-d}^{j}a_{n}e^{n\eta }}{\sum_{s=-q}^{p}b_{m}e^{m\eta }% }=\frac{a_{j}e^{j\eta }+\ldots +a_{-d}e^{-d\eta }}{b_{p}e^{p\eta }+\ldots +a_{-q}e^{-q\eta }} \end{array}$$(10)

where p, q, j and d are positive integers that can be examined later, a_n and b_m are unrecognized constants. To calculate the values of j and p, highest order the linear term of Eq. (9) is equalized with the highest order nonlinear term. By using the same procedure, the values for q and d, can be calculated by balancing the lowest order linear term of Eq. (9) with lowest order nonlinear term. As a result we can acquire the traveling wave solutions of the considered Eq. (6)

Formerly, a perturbation based algorithm has been introduced by Aksoy and Pakdemirli [3]. In the method, an iterative algorithm is proposed using the perturbation expansion. Previously, this method is implemented on ordinary FDEs [27], fractional-integro differential equations [28] and systems of FDEs [29].

In this article, the most basic PIA, PIA(1, 1) is used to attain approximate solutions of FPDEs. For this purpose, one we consider the correction term in the perturbation expansion and correction terms of first derivatives in the Taylor series expansion [3, 4].

To describe the main idea of PIA, take the FPDE

Fuα<!-- a -->(x,t),u(x,t),ux(x,t),uxx(x,t),…<!-- ? -->,ε<!-- e -->=0, $$\begin{array}{} \displaystyle F\left( u^{\alpha }(x,t),u(x,t),u_{x}(x,t),u_{xx}(x,t),\dots ,\varepsilon \right) =0, \end{array}$$(11)

where ε is assumed as an artificially small parameter. The perturbation expansions with one correction terms are

un+1=un+ε<!-- e -->ucnTα<!-- a -->(un+1)=Tα<!-- a -->(un)+ε<!-- e -->uc′<!-- ' -->n $$\begin{array}{} \begin{split} \displaystyle u_{n+1}& =u_{n}+\varepsilon {\left( u_{c}\right) }_{n} \\ T_{\alpha }(u_{n+1})& =T_{\alpha }(u_{n})+\varepsilon {\left( u_{c}^{\prime }\right) }_{n} \end{split} \end{array}$$(12)

Subrogating (12) into (11) and expanding in the Taylor series form for only first order derivatives yields

Funα<!-- a -->,un,(un)x,(un)xx,…<!-- ? -->,0+Fuunα<!-- a -->,un,(un)x,(un)xx,…<!-- ? -->,0ε<!-- e -->ucn+Fuα<!-- a -->unα<!-- a -->,un,(un)x,(un)xx,…<!-- ? -->,0ε<!-- e -->uc(α<!-- a -->)n+Fuxunα<!-- a -->,un,(un)x,(un)xx,…<!-- ? -->,0ε<!-- e -->(uc)xn+Fuxxunα<!-- a -->,un,(un)x,(un)xx,…<!-- ? -->,0ε<!-- e -->(uc)xxn+⋯<!-- ? -->+Fε<!-- e -->unα<!-- a -->,un,unx,unxx,…<!-- ? -->,0ε<!-- e -->=0 $$\begin{array}{} \begin{split} \displaystyle F\left( u_{n}^{\left( \alpha \right) },u_{n},{(u_{n})}_{x},{(u_{n})}% _{xx},\dots ,0\right) & +F_{u}\left( u_{n}^{\left( \alpha \right) },u_{n},{% (u_{n})}_{x},{(u_{n})}_{xx},\dots ,0\right) \varepsilon {\left( u_{c}\right) }_{n} \\ & +F_{u^{\left( \alpha \right) }}\left( u_{n}^{\left( \alpha \right) },u_{n},% {(u_{n})}_{x},{(u_{n})}_{xx},\dots ,0\right) \varepsilon {\left( u_{c}^{(\alpha )}\right) }_{n} \\ & +F_{u_{x}}\left( u_{n}^{\left( \alpha \right) },u_{n},{(u_{n})}_{x},{% (u_{n})}_{xx},\dots ,0\right) \varepsilon {\left( {(u_{c})}_{x}\right) }_{n} \\ & +F_{u_{xx}}\left( u_{n}^{\left( \alpha \right) },u_{n},{(u_{n})}_{x},{% (u_{n})}_{xx},\dots ,0\right) \varepsilon {\left( {(u_{c})}_{xx}\right) }% _{n}+\cdots \\ & +F_{\varepsilon }\left( u_{n}^{\left( \alpha \right) },u_{n},{\left( u_{n}\right) }_{x},{\left( u_{n}\right) }_{xx},\dots ,0\right) \varepsilon =0 \end{split} \end{array}$$(13)

or

uc(α<!-- a -->)n∂<!-- ? -->F∂<!-- ? -->u(α<!-- a -->)+ucn∂<!-- ? -->F∂<!-- ? -->u+(uc)xn∂<!-- ? -->F∂<!-- ? -->ux+(uc)xxn∂<!-- ? -->F∂<!-- ? -->uxx+⋯<!-- ? -->+∂<!-- ? -->F∂<!-- ? -->ε<!-- e -->+Fε<!-- e -->=0. $$\begin{array}{} \begin{split} \displaystyle {\left( u_{c}^{(\alpha )}\right) }_{n}\frac{\partial F}{\partial u^{(\alpha )}}+{\left( u_{c}\right) }_{n}\frac{\partial F}{\partial u}+{\left( {(u_{c})}% _{x}\right) }_{n}\frac{\partial F}{\partial u_{x}}+{\left( {(u_{c})}% _{xx}\right) }_{n}\frac{\partial F}{\partial u_{xx}}+\dots +\frac{\partial F% }{\partial \varepsilon }+\frac{F}{\varepsilon }=0. \end{split} \end{array}$$(14)

Rewriting (14) gives the subsequent PIA(1, 1) iteration formula

ut(x,t)+FuFutu(x,t)=−<!-- - -->Fε<!-- e -->+Fε<!-- e -->Fut. $$\begin{array}{} \begin{split} \displaystyle u_{t}(x,t)+\frac{F_{u}}{F_{u_{t}}}u(x,t)=-\frac{F_{\varepsilon }+\frac{F}{% \varepsilon }}{F_{u_{t}}}. \end{split} \end{array}$$(15)

In this expansion, all of the derivatives are evaluated at ε = 0. Using an initial function u₀(x, t), firstly the correction term (u_c)₀(x, t) is computed. Subrogating it into (12) gives the first approximate result u₁(x, t). Similar procedure is applied until obtaining the other approximations.

Think of the fractional coupled Burgers’ equation [24] as

Dtα<!-- a -->u+uDxu+vDyu−<!-- - -->1ℜ<!-- R -->(Dx2u+Dy2u)=0,Dtα<!-- a -->v+uDxv+vDyv−<!-- - -->1ℜ<!-- R -->(Dx2v+Dy2v)=0, $$\begin{array}{} \begin{split} \displaystyle D_{t}^{\alpha }{u}+uD_{x}{u}+vD_{y}{u}-\frac{1}{\Re }(D_{x}^{2}{u}+D_{y}^{2}{% u})& ={0,} \\ D_{t}^{\alpha }{v}+uD_{x}{v}+vD_{y}{v}-\frac{1}{\Re }(D_{x}^{2}{v}+D_{y}^{2}{% v})& =0,\, \end{split} \end{array}$$(16)

where 0 < α ≤ 1, 𝔎, is Reynolds number, u = u(x, y, t) and v = v(x, y, t). By the help of the chain rule [1] and the wave transform η<!-- ? -->=kx+wy+ctα<!-- a -->α<!-- a -->, $\begin{array}{} \eta =kx+wy+c\frac{t^{\alpha }}{\alpha }, \end{array}$ we obtain

cU′<!-- ' -->(η<!-- ? -->)+ku(η<!-- ? -->)U′<!-- ' -->(η<!-- ? -->)+wV(η<!-- ? -->)U′<!-- ' -->(η<!-- ? -->)−<!-- - -->1ℜ<!-- R -->(k2+w2)U′<!-- ' -->′<!-- ' -->(η<!-- ? -->)=0cV′<!-- ' -->(η<!-- ? -->)+kuV′<!-- ' -->(η<!-- ? -->)+wvV′<!-- ' -->(η<!-- ? -->)−<!-- - -->1ℜ<!-- R -->(k2+w2)V′<!-- ' -->′<!-- ' -->(η<!-- ? -->)=0. $$\begin{array}{} \begin{split} \displaystyle & cU^{\prime }(\eta )+ku(\eta )U^{\prime }(\eta )+wV(\eta )U^{\prime }(\eta )-\frac{1}{\Re }(k^{2}+w^{2})U^{\prime \prime }(\eta )=0 \\ & cV^{\prime }(\eta )+kuV^{\prime }(\eta )+wvV^{\prime }(\eta )-\frac{1}{\Re }(k^{2}+w^{2})V^{\prime \prime }(\eta )=0. \end{split} \end{array}$$(17)

Now assume that the solution of (17) can be described as

U(η<!-- ? -->)=acecη<!-- ? -->+…<!-- ? -->+a−<!-- - -->de−<!-- - -->dη<!-- ? -->bpepη<!-- ? -->+…<!-- ? -->+b−<!-- - -->qe−<!-- - -->qη<!-- ? -->. $$\begin{array}{} \begin{split} \displaystyle U(\eta )=\frac{a_{c}e^{c\eta }+\ldots +a_{-d}e^{-d\eta }}{b_{p}e^{p\eta }+\ldots +b_{-q}e^{-q\eta }}. \end{split} \end{array}$$(18)

V(η<!-- ? -->)=dsesη<!-- ? -->+…<!-- ? -->+d−<!-- - -->ne−<!-- - -->nη<!-- ? -->flelη<!-- ? -->+…<!-- ? -->+f−<!-- - -->re−<!-- - -->rη<!-- ? -->. $$\begin{array}{} \begin{split} \displaystyle V(\eta )=\frac{d_{s}e^{s\eta }+\ldots +d_{-n}e^{-n\eta }}{f_{l}e^{l\eta }+\ldots +f_{-r}e^{-r\eta }}. \end{split} \end{array}$$(19)

Using (18), (19) and (17) led to c = p, d = q, s = l and n = r. For convenience lets assume all the coefficients c = p = s = l = n = r = 1. Now rewriting u(η) and v(η) due to above assumptions

U=a1eη<!-- ? -->+a0+a−<!-- - -->1e−<!-- - -->η<!-- ? -->b1eη<!-- ? -->+b0+b−<!-- - -->1e−<!-- - -->η<!-- ? -->V=d1eη<!-- ? -->+d0+d−<!-- - -->1e−<!-- - -->η<!-- ? -->f1eη<!-- ? -->+f0+f−<!-- - -->1e−<!-- - -->η<!-- ? -->. $$\begin{array}{} \begin{split} \displaystyle & U=\frac{a_{{1}}{e^{\eta }}+a_{{0}}+a_{{-1}}{e^{-\eta }}}{b_{{1}}{e^{\eta }}% +b_{{0}}+b_{{-1}}{e^{-\eta }}} \\ & V={\frac{d_{{1}}{e^{\eta }}+d_{{0}}+d_{{-1}}{e^{-\eta }}}{f_{{1}}{e^{\eta }% }+f_{{0}}+f_{{-1}}{e^{-\eta }}}}. \end{split} \end{array}$$(20)

Substituting the equations (20) into (17) and equalizing the coefficients of e^nη yields a system of algebraic equations. Solving the system with respect to the constants expressed above we can handle the following solutions

c=−<!-- - -->ℜ<!-- R -->wd0b0+k2f0b0+w2f0b0+ka0f0ℜ<!-- R -->b0f0,a−<!-- - -->1=f−<!-- - -->1ℜ<!-- R -->wd0b0+2k2b0f−<!-- - -->1f0+2w2b0f−<!-- - -->1f0+ℜ<!-- R -->ka0f−<!-- - -->1f0−<!-- - -->ℜ<!-- R -->b0wd−<!-- - -->1f0ℜ<!-- R -->kf02,a1=0,b−<!-- - -->1=b0f−<!-- - -->1f0,a1=0,b1=0,d1=0,f1=0. $$\begin{array}{} \begin{split} \displaystyle & c=-{\frac{\Re wd_{{0}}b_{{0}}+{k}^{2}f_{{0}}b_{{0}}+{w}^{2}f_{{0}}b_{{0}% }+ka_{{0}}f_{{0}}}{\Re b_{{0}}f_{{0}}}}, \\ & a_{-1}={\frac{f_{{-1}}\Re wd_{{0}}b_{{0}}+2\,{k}^{2}b_{{0}}f_{{-1}}f_{{0}% }+2\,{\ w}^{2}b_{{0}}f_{{-1}}f_{{0}}+\Re ka_{{0}}f_{{-1}}f_{{0}}-\Re b_{{0}% }wd_{{-1}}f_{{0}}}{\Re k{f_{{0}}}^{2}}}, \\ & a_{1}=0,b_{-1}={\frac{b_{{0}}f_{{-1}}}{f_{{0}}}}, \\ & a_{1}=0,b_{1}=0,d_{1}=0,f_{1}=0. \end{split} \end{array}$$(21)

So the solutions can be obtained as

u(x,y,t)=a0+AeBb0+b0f−<!-- - -->1eBf0−<!-- - -->1ℜ<!-- R -->kf02 $$\begin{array}{} \begin{split} \displaystyle u(x,y,t)=\frac{a_{{0}}+A {e^{B}}b_{0}+b_{0}f_{{-1}}{e^{B}}{f_{{0}}}^{-1}}{% \Re {k}{f_{{0}}}^{2}} \end{split} \end{array}$$(22)

and

v(x,y,t)=d0+d−<!-- - -->1eBf0+f−<!-- - -->1eB $$\begin{array}{} \begin{split} \displaystyle v(x,y,t)=\frac{d_{{0}}+d_{{-1}}{e^{B}}}{f_{{0}}+f_{{-1}}{e^{B}}} \end{split} \end{array}$$(23)

where

A=f−<!-- - -->1ℜ<!-- R -->wd0b0+2k2b0f−<!-- - -->1f0+2w2b0f−<!-- - -->1f0+ℜ<!-- R -->ka0f−<!-- - -->1f0−<!-- - -->ℜ<!-- R -->b0wd−<!-- - -->1f0B=−<!-- - -->kx−<!-- - -->wy+ℜ<!-- R -->wd0b0+k2f0b0+w2f0b0+ℜ<!-- R -->ka0f0tα<!-- a -->ℜ<!-- R -->b0f0α<!-- a --> $$\begin{array}{} \begin{split} \displaystyle & A=\left( f_{{-1}}\Re wd_{{0}}b_{{0}}+2{k}^{2}b_{{0}}f_{{-1}}f_{{0}% }+2\,{w}^{2}b_{{0}}f_{{-1}}f_{{0}}+\Re ka_{{0}}f_{{-1}}f_{{0}}-\Re b_{{0}% }wd_{{-1}}f_{{0}}\right) \\ & B=-kx-wy+{\frac{\left( \Re wd_{{0}}b_{{0}}+{k}^{2}f_{{0}}b_{{0}}+{w}^{2}f_{% {0}}b_{{0}}+\Re ka_{{0}}f_{{0}}\right) {t}^{\alpha }}{\Re b_{{0}}f_{{0}% }\alpha }} \end{split} \end{array}$$(24)

Regard the system (16) with the the conditions u(x,y,0)=1+14ℜ<!-- R -->ℜ<!-- R -->+8e−<!-- - -->x−<!-- - -->y1+12e−<!-- - -->x−<!-- - -->yandv(x,y,0)=1+e−<!-- - -->x−<!-- - -->y2+e−<!-- - -->x−<!-- - -->y. $\begin{array}{} {u}(x,y,0)=\frac{% 1+\frac{1}{4\Re }\left( \Re +8\right) {e^{-x-y}}}{1+\frac{1}{2}{e^{-x-y}}}~~ \text{and}~~ {v}(x,y,0)=\frac{1+{e^{-x-y}}}{2+{e^{-x-y}}}. \end{array}$ For the values k = 1, w = 1, f₀ = 2, d₀ = 1, b₀ = 1, d_–1 = 1, b_–1 = 1, f_–1 = 1 and a₀ = 1, we can acquire the exact solutions as

u(x,y,t)=1+14ℜ<!-- R -->+8e−<!-- - -->x−<!-- - -->y+3ℜ<!-- R -->+4tα<!-- a -->2ℜ<!-- R -->α<!-- a -->1+12e−<!-- - -->x−<!-- - -->y+3ℜ<!-- R -->+4tα<!-- a -->2Rα<!-- a -->,v(x,y,t)=1+e−<!-- - -->x−<!-- - -->y+3ℜ<!-- R -->+4tα<!-- a -->2ℜ<!-- R -->α<!-- a -->2+e−<!-- - -->x−<!-- - -->y+3ℜ<!-- R -->+4tα<!-- a -->2ℜ<!-- R -->α<!-- a -->. $$\begin{array}{} \begin{split} \displaystyle u(x,y,t)& =\frac{1+\frac{1}{4}\,\left( \Re +8\right) {e^{-x-y+\,{\frac{% \left( 3\,\Re +4\right) {t}^{\alpha }}{2\Re \alpha }}}}}{1+\frac{1}{2}\,{% e^{-x-y+\,{\frac{\left( 3\,\Re +4\right) {t}^{\alpha }}{2\mathfrak{R}\alpha }% }}}}, \\ v(x,y,t)& =\frac{1+{e^{-x-y+\,{\frac{\left( 3\,\Re +4\right) {t}^{\alpha }}{% 2\Re \alpha }}}}}{2+{e^{-x-y+\,{\frac{\left( 3\,\Re +4\right) {t}^{\alpha }}{% 2\Re \alpha }}}}}. \end{split} \end{array}$$(25)

Now we introduce a small perturbation parameter ε to the system and rewrite the equations as

Dtα<!-- a -->u+ε<!-- e -->uDxu+ε<!-- e -->vDyu−<!-- - -->1ℜ<!-- R -->(Dx2u+Dy2u)=0,Dtα<!-- a -->v+ε<!-- e -->uDxv+ε<!-- e -->vDyv−<!-- - -->1ℜ<!-- R -->(Dx2v+Dy2v)=0, $$\begin{array}{} \begin{split} \displaystyle D_{t}^{\alpha }{u}+{\varepsilon }uD_{x}{u}+{\varepsilon }vD_{y}{u}-\frac{1}{% \Re }(D_{x}^{2}{u}+D_{y}^{2}{u})& ={0,} \\ D_{t}^{\alpha }{v}+{\varepsilon }uD_{x}{v}+{\varepsilon }vD_{y}{v}-\frac{1}{% \Re }(D_{x}^{2}{v}+D_{y}^{2}{v})& =0,\, \end{split} \end{array}$$(26)

Therefore, terms in formula (15)turn into

F=t1−<!-- - -->α<!-- a -->(un)t(x,y,t)−<!-- - -->1ℜ<!-- R -->(un)xx(x,y,t)+(un)yy(x,y,t),Fu=0,Fut=t1−<!-- - -->α<!-- a -->,Fε<!-- e -->=un(x,y,t)unx(x,y,t)+vn(x,y,t)(un)y(x,y,t) $$\begin{array}{} \begin{split} \displaystyle F& =t^{1-\alpha }(u_{n})_{t}(x,y,t)-\frac{1}{\Re }\left( (u_{n})_{xx}(x,y,t)+(u_{n})_{yy}(x,y,t)\right) ,\text{ }F_{{u}}=0,\quad \\ F_{{u}_{t}}& =t^{1-\alpha },\text{ }F_{\varepsilon }=u_{n}(x,y,t)\left( {u}% _{n}\right) _{x}(x,y,t)+v_{n}(x,y,t)(u_{n})_{y}(x,y,t) \end{split} \end{array}$$(27)

and

F=t1−<!-- - -->α<!-- a -->(vn)t(x,y,t)−<!-- - -->1ℜ<!-- R -->(vn)xx(x,y,t)+(vn)yy(x,y,t),Fv=0,Fvt=t1−<!-- - -->α<!-- a -->,Fε<!-- e -->=un(x,t)(vn)x(x,t)+vn(x,t)(vn)y(x,t) $$\begin{array}{} \begin{split} \displaystyle F& =t^{1-\alpha }(v_{n})_{t}(x,y,t)-\frac{1}{\Re }\left( (v_{n})_{xx}(x,y,t)+(v_{n})_{yy}(x,y,t)\right) ,\text{ }F_{{v}}=0, \\ F_{v_{t}}& =t^{1-\alpha },\text{ }F_{\varepsilon }=u_{n}(x,t)(v_{n})_{x}(x,t)+v_{n}(x,t)(v_{n})_{y}(x,t) \end{split} \end{array}$$(28)

Subrogating above terms in the iteration formula (15) gives the subsequent partial differential equations

ε<!-- e -->ℜ<!-- R -->tuct=−<!-- - -->ℜ<!-- R -->tunt+tα<!-- a -->((un)xx−<!-- - -->ε<!-- e -->ℜ<!-- R -->(un(un)x+vn(un)y)+(un)yy) $$\begin{array}{} \begin{split} \displaystyle \varepsilon \Re \mathfrak{t}\left( u_{c}\right) _{t}=-\Re t\left( u_{n}\right) _{t}+t^{\alpha }((u_{n})_{xx}-\varepsilon \Re (u_{n}(u_{n})_{x}+v_{n}(u_{n})_{y})+(u_{n})_{yy}) \end{split} \end{array}$$(29)

and

ε<!-- e -->ℜ<!-- R -->tvct=−<!-- - -->ℜ<!-- R -->tvnt+tα<!-- a -->((vn)xx−<!-- - -->ε<!-- e -->ℜ<!-- R -->(un(vn)x+vn(vn)y)+(vn)yy) $$\begin{array}{} \begin{split} \displaystyle \varepsilon \Re \mathfrak{t}\left( v_{c}\right) _{t}=-\Re t\left( v_{n}\right) _{t}+t^{\alpha }((v_{n})_{xx}-\varepsilon \Re (u_{n}(v_{n})_{x}+v_{n}(v_{n})_{y})+(v_{n})_{yy}) \end{split} \end{array}$$(30)

Beginning with the initial functions

u(x,y,0)=1+14ℜ<!-- R -->ℜ<!-- R -->+8e−<!-- - -->x−<!-- - -->y1+12e−<!-- - -->x−<!-- - -->yandv(x,y,0)=1+e−<!-- - -->x−<!-- - -->y2+e−<!-- - -->x−<!-- - -->y $$\begin{array}{} \begin{split} \displaystyle {u}(x,y,0)=\frac{1+\frac{1}{4\Re }\left( \Re +8\right) {e^{-x-y}}}{1+\frac{1% }{2}{e^{-x-y}}}~\text{and}~{v}(x,y,0)=\frac{1+{e^{-x-y}}}{2+{e^{-x-y}}} \end{split} \end{array}$$(31)

and using (15), the numerical results are obtained for n = 0, 1, 2, … respectively.

u1(x,y,t)=4ℜ<!-- R -->ex+y+ℜ<!-- R -->+82ℜ<!-- R -->2ex+y+1−<!-- - -->(ℜ<!-- R -->−<!-- - -->8)(3ℜ<!-- R -->+4)tα<!-- a -->ex+y2α<!-- a -->ℜ<!-- R -->22ex+y+12 $$\begin{array}{} \begin{split} \displaystyle u_{1}(x,y,t) &=&\frac{4\Re e^{x+y}+\Re +8}{2\Re \left( 2e^{x+y}+1\right) }-% \frac{(\Re -8)(3\Re +4)t^{\alpha }e^{x+y}}{2\alpha \Re ^{2}\left( 2e^{x+y}+1\right) ^{2}} \end{split} \end{array}$$(32)

v1(x,y,t)=(3ℜ<!-- R -->+4)tα<!-- a -->ex+y2α<!-- a -->ℜ<!-- R -->2ex+y+12+ex+y+12ex+y+1 $$\begin{array}{} \begin{split} \displaystyle v_{1}(x,y,t) &=&\frac{(3\Re +4)t^{\alpha }e^{x+y}}{2\alpha \Re \left( 2e^{x+y}+1\right) ^{2}}+\frac{e^{x+y}+1}{2e^{x+y}+1} \end{split} \end{array}$$(33)

u2(x,y,t)=−<!-- - -->(ℜ<!-- R -->−<!-- - -->8)(3ℜ<!-- R -->+4)2t2α<!-- a -->ex+y2ex+y−<!-- - -->13α<!-- a -->ℜ<!-- R -->2ex+y+12+16tα<!-- a -->ex+y24α<!-- a -->3ℜ<!-- R -->42ex+y+15−<!-- - -->(ℜ<!-- R -->−<!-- - -->8)(3ℜ<!-- R -->+4)tα<!-- a -->ex+y2α<!-- a -->2ℜ<!-- R -->ex+y+ℜ<!-- R -->2+8−<!-- - -->ℜ<!-- R -->2ℜ<!-- R -->2ex+y+1+1 $$\begin{array}{} \begin{split} \displaystyle u_{2}(x,y,t) &=&-\frac{(\Re -8)(3\Re +4)^{2}t^{2\alpha }e^{x+y}\left( 2e^{x+y}-1\right) \left( 3\alpha \Re \left( 2e^{x+y}+1\right) ^{2}+16t^{\alpha }e^{x+y}\right) }{24\alpha ^{3}\Re ^{4}\left( 2e^{x+y}+1\right) ^{5}} \notag \\ &&-\frac{(\Re -8)(3\Re +4)t^{\alpha }e^{x+y}}{2\alpha \left( 2\Re e^{x+y}+\Re \right) ^{2}}+\frac{8-\Re }{2\Re \left( 2e^{x+y}+1\right) }+1 \end{split} \end{array}$$(34)

v2(x,y,t)=2(3ℜ<!-- R -->+4)2t3α<!-- a -->e2x+2y2ex+y−<!-- - -->13α<!-- a -->3ℜ<!-- R -->32ex+y+15+(3ℜ<!-- R -->+4)2t2α<!-- a -->ex+y2ex+y−<!-- - -->18α<!-- a -->2ℜ<!-- R -->22ex+y+13+(3ℜ<!-- R -->+4)tα<!-- a -->ex+y2α<!-- a -->ℜ<!-- R -->2ex+y+12+ex+y+12ex+y+1 $$\begin{array}{} \begin{split} \displaystyle v_{2}(x,y,t) &=&\frac{2(3\Re +4)^{2}t^{3\alpha }e^{2x+2y}\left( 2e^{x+y}-1\right) }{3\alpha ^{3}\Re ^{3}\left( 2e^{x+y}+1\right) ^{5}}+\frac{% (3\Re +4)^{2}t^{2\alpha }e^{x+y}\left( 2e^{x+y}-1\right) }{8\alpha ^{2}\Re ^{2}\left( 2e^{x+y}+1\right) ^{3}} \notag \\ &&+\frac{(3\Re +4)t^{\alpha }e^{x+y}}{2\alpha \Re \left( 2e^{x+y}+1\right) ^{2}}+\frac{e^{x+y}+1}{2e^{x+y}+1} \end{split} \end{array}$$(35)

Similarly, the fourth order solutions u₄(x, y, t) and v₄(x, y, t) are calculated. In Table 1 and 2, the fourth order PIA numerical approximate solutions are compared to exact solutions. Also the absolute errors are calculated for changing values of α and x. The results indicate the reliability of PIA. Besides using Figures 1 – 6, figures for the solutions of PIA are illustrated for changing α values. They exhibit that PIA produces highly approximate results. It is also obvious that further iterations would generate convenient solutions.

Fig. 1

Fig. 2

Fig. 3

Fig. 4

Fig. 5

Fig. 6

Table 1

PIA (u₄(x, y, t)) and exact solution values with absolute errors for y = 1, t = 0.1 and ℜ = 100.

Table 2

PIA (v₄(x, y, t)) and exact solution values with absolute errors for y = 1, t = 0.1 and ℜ = 100.