# Approximate Solutions to Time-Fractional Schrödinger Equation via Homotopy Analysis Method

Authors: Najeeb Alam Khan; Muhammad Jamil; Asmat Ara
Journal: ISRN Mathematical Physics
 (2012)

---

## Abstract

We construct the approximate solutions of the time-fractional Schrödinger equations, with zero and nonzero trapping potential, by homotopy analysis method (HAM). The fractional derivatives, in the Caputo sense, are used. The method is capable of reducing the size of calculations and handles nonlinear-coupled equations in a direct manner. The results show that HAM is more promising, convenient, efficient and less computational than differential transform method (DTM), and easy to apply in spaces of higher dimensions as well.


---

## Body

## 1. Introduction

The theory of derivatives of fractional (nonintegers) orders stimulates considerable interest in the areas of mathematics, physics, engineering, and other sciences. Fractional derivatives [1–5] provide an excellent tool for the description of memory and hereditary properties of various material and processes. The beauty of this subject is that a fractional derivative is not a local point property. This considers the history and nonlocal distributed effects. Perhaps, this subject translates the reality of nature better. Application of fractional calculus are found in different areas of sciences such as physics, continuum mechanics, signal processing, electromagnetics, and bioengineering. The electrical properties of nerve cell membranes and the propagation of electrical signals are well characterized by differential equations of fractional order. The fractional differential equations (FDE) [6–11] appear more and more frequently in different research areas and engineering applications. Exact solution of nonlinear partial differential equations has become one of the central themes of perpetual interest in mathematical physics. In order to better understand these phenomena as well as further apply them in the practical life, it is important to seek their exact solutions.A great deal of effort has also been expanded in attempting to find robust and stable numerical and analytical methods for solving fractional differential equations of physical interest. These methods include Laplace transform method, Fourier transform method, finite difference method [12], fractional linear multistep methods, Adomian decomposition method (ADM) [13], variational iteration method (VIM) [14], homotopy perturbation method (HPM) [15], differential transform method [16], and homotopy analysis method (HAM) [17–20].The aim of this paper is to solve the Schrödinger equation with fractional order using the homotopy analysis method. By introducing an embedding parameterq the nonlinear fractional differential equation is converted to a linear fractional differential equation at q=0. When q evolves, the differential equation becomes the original one at q=1. The method has been used in a variety of problems and the details can be found in Liao’s book [17]. The method gives rapidly convergent successive approximations of the exact solution if such solution exists; otherwise, a few approximations can be used for numerical purposes. Wang [12] presented the numerical solution Schrödinger equations by means of finite difference scheme. Khuri [13] applied ADM to obtain the solution of cubic Schrödinger equations. Wazwaz [14] presented the exact solution of the linear and nonlinear one-dimensional Schrödinger equations by VIM. Recently, Ravi Kanth [16] and his colleague presented the exact solution of the linear and nonlinear Schrödinger equations by differential transformation method (DTM).The aim of this paper is to investigate the approximate solutions of the time-fractional Schrödinger equations, with zero and non-zero trapping potential, by means of HAM. The convergent region is then obtained by looking at the real and imaginary parts of the series in plot.


## 2. Fractional Schrödinger Equations and Preliminaries

The time-fractional Schrödinger equations (FSE) has the following form:(2.1)i∂αψ(X,t)∂tα+12∇2ψ(X,t)-Vd(X)ψ-βd|ψ|2ψ=0;t≥0,ψ(X,0)=ψ0(X),X∈Rd,
where Vd(X) is the trapping potential and βd is a real constant. The physical model (2.1) and its generalized forms occur in various areas of physics, including nonlinear optics, plasma physics, superconductivity, and quantum mechanics.We give some basic definitions, notations, and properties of the fractional calculus theory, which will be used later in this work.Definition 2.1.
The Riemann-Liouville fractional integral operator of orderμ on the usual Lebesgue space L1[a,b] is given by
(2.2)Jμf(x)=1Γ(μ)∫0x(x-t)μ-1f(t)dt;μ>0,(2.3)J0f(x)=f(x).
This integral operator has the following properties:(i)JαJβ=Jα+β=JβJα,α,β>0,(ii)Jα(x-a)γ=Γ(γ+1)Γ(α+γ+1)(x-a)α+γ,α>0,γ>-1.Definition 2.2.
The Caputo definition of fractal derivative operator is given by(2.4)Dμf(x)=1Γ(m-μ)∫0t(x-τ)m-μ-1f(m)(τ)dτm-1<μ≤m,m∈N,x>0.
It has the following two basic properties:
(2.5)DμJμf(x)=f(x),
and
(2.6)JμDμf(x)=f(x)-∑k=0m-1f(k)(0+)(x)kk!,x>0.
The Caputo fractional derivative is considered here because it allows traditional initial and boundary conditions to be included in the formulation of the problem.


## 3. Basic Idea of HAM

In this paper, we apply the HAM to the linear and nonlinear problems to be discussed. In order to show the basic idea of HAM, consider the following nonlinear fractional differential equation:(3.1)NF[ψ(X,t;q)]=0,where𝒩ℱ is a nonlinear fractional operator, X and t denote the independent variables, and ψ is an unknown function. By means of the HAM, we first construct the so-called zeroth-order deformation equation(3.2)(1-q)£F[ϕ(X,t;q)-ψ0(X,t)]=qħNF[ϕ(X,t;q)],whereq∈[0,1] is the embedding parameter, ħ≠0 is an auxiliary parameter, £ℱ is an auxiliary linear operator, ϕ(X,t;q)  are unknown functions, and ψ0(X,t) are initial guesses of ϕ(X,t;q). It is obvious that for q=0 and q=1, (3.3) becomes(3.3)ϕ(X,t;0)=ψ0(X,t),ϕ(X,t;1)=ψ(X,t),respectively. Thus asq=0 increases from 0 to 1, the solution ϕ(X,t;q) varies from the initial guess ψ0(X,t) to the exact solution ψ(X,t). Expanding ϕ(X,t;q) in Taylor series with respect to q, one has(3.4)ϕ(X,t;q)=ψ0(X,t)+∑m=1+∞ψm(X,t)qm,where(3.5)ψm(X,t)=1m!∂mϕ(X,t;q)∂qm|q=0.The convergence of the series (3.4) depends upon the auxiliary parameter ħ. If it is convergent at q=1, one has(3.6)ψ(X,t)=ψ0(X,t)+∑m=1+∞ψm(X,t),which must be one of the solutions of the original nonlinear equations, as proved by Liao [17]. Define the vectors(3.7)ψn⇀={ψ0(X,t),ψ1(X,t),ψ2(X,t),…,ψn(X,t)}.
Differentiating the zeroth-order deformation (3.2) m-times with respect to q and then dividing them by m!,  and finally setting q=0, we get the following mth-order deformation equation:(3.8)£F[ψm(X,t)-χmψm-1(X,t)]=ħRm(ψ⃗m-1),where(3.9)Rm(ψ⃗m-1)=1(m-1)!∂m-1NF[ϕ(X,t;q)]∂qm-1|q=0,χm={0;m=1,1;m>1.It should be emphasized thatψm(X,t;q) for m≥1 is governed by the linear equations (3.8) with boundary conditions that come from the original problem, which can be solved by the symbolic computation software MATHEMATICA. The success of the technique is based on the proper selection of the initial guess.


## 4. Implementation of the Method

In this section, we introduce the above reliable approach, in a realistic and efficient way, to handle nonlinear Schrödinger equation with time-fractional derivative.


### 4.1. Nonlinear-Time Fractional Schrödinger Equation (NLFSE)

Settingψ(X,t)=v(X,t)+iw(X,t) and ψ(X,0)=v(X,0)+iw(X,0) in (2.1) leads to the following coupled system of equations(4.1)Dαv+12∇2w-Vdw-βd(wv2+w3)=0,(4.2)Dαw-12∇2v-Vdv+βd(v3+w2v)=0,subject to the initial conditions(4.3)v(X,0)=v0=F(X),w(X,0)=w0=G(X).Equations (4.1) and (4.3) suggest that we define the nonlinear operator as(4.4)NF1[φ1(X,t;q),φ2(X,t;q)]=Dαφ1(X,t;q)+12∇2φ2(X,t;q)-Vdφ2(X,t;q)-φ2(X,t;q)φ12(X,t;q)-φ23(X,t;q),(4.5)NF2[φ1(X,t;q),φ2(X,t;q)]=Dαφ2(X,t;q)-12∇2φ1(X,t;q)+Vdφ1(X,t;q)+φ13(X,t;q)+φ22(X,t;q)φ1(X,t;q),and the linear operator(4.6)£F[ϕ(X,t;q)]=Dα[ϕ(X,t;q)]with the property£F[c1(X)]=0£F[vm(X,t)-χmvm-1(X,t)]=ħR1,m[v⃗m-1,w⃗m-1](4.7)£F[wm(X,t)-χmwm-1(X,t)]=ħR2,m[v⃗m-1,w⃗m-1],where(4.8)R1,m[v⃗m-1,w⃗m-1]=Dαvm-1+12∇2wm-1-Vdwm-1-βd(∑k=0m-1vm-k-1∑l=0kvk-lwl+∑k=0m-1wm-k-1∑l=0kwk-lwl)R2,m[v⃗m-1,w⃗m-1]=Dαwm-1-12∇2vm-1+Vdvm-1+βd(∑k=0m-1wm-k-1∑l=0kwk-lvl+∑k=0m-1vm-k-1∑l=0kvk-lvl).Obviously, the solution of themth-order deformation (4.8) for m ≥ 1 becomes(4.9)vm(X,t)=χmvm-1(X,t)+ħJαR1,m[v⃗m-1,w⃗m-1]wm(X,t)=χmwm-1(X,t)+ħJαR2,m[v⃗m-1,w⃗m-1].Example 4.1.
We first consider the one-dimensional NLFSE with zero trapping potential (i.e.,Vd(x)=0) and βd=-1(4.10)iDtαψ(x,t)+12∂2ψ(x,t)∂x2+|ψ|2ψ=0,t≥0,0<α≤1.
subject to the initial condition ψ(x,0)=eix.
Solving the above equations(4.11)v0=cos⁡x,w0=sin⁡xv1=ħtαsin⁡x2Γ(α+1),w2=-ħtαcos⁡x2Γ(α+1)v2=ħ(ħ+1)tαsin⁡x2Γ(α+1)-ħ2t2αcos⁡x4Γ(2α+1),w2=-ħ(ħ+1)tαcos⁡x2Γ(α+1)-ħ2t2αsin⁡x4Γ(2α+1)v3=(ħ(ħ+1)2tα2Γ(α+1)-ħ3t3α(5(Γ(α+1))2-2Γ(2α+1))8(Γ(α+1))2Γ(3α+1))sin⁡x-ħ2(ħ+1)t2αcos⁡x2Γ(2α+1),w3=(-ħ(ħ+1)2tα2Γ(α+1)+ħ3t3α(5(Γ(α+1))2-2Γ(2α+1))8(Γ(α+1))2Γ(3α+1))cos⁡x-ħ2(ħ+1)t2αsin⁡x2Γ(2α+1),and so on, in this manner the rest of the components can be obtained. Therefore, the approximate solution is
(4.12)ψ(x,t)=∑n=07(vn+iwn).
The exact solution of (4.10) for α=1 is ψ=ei(x+(t/2)). When ħ=-1, α=1, the solution obtained by [13–16] is recovered as a special case.Example 4.2.
Consider the one-dimensional NLFSE with trapping potential, that is,Vd(x)=cos2x and βd=-1(4.13)iDtαψ(x,t)+12∂2ψ(x,t)∂x2-ψcos⁡2⁡x-|ψ|2ψ=0,t≥0,0<α≤1.
subject to the initial condition ψ(x,0)=sinx.
Solving above equations we obtain the solution in a series form(4.14)v0=sin⁡x,w0=0v1=0,w1=3ħtαsin⁡x2Γ(α+1)v2=-9ħ2t2αsin⁡x4Γ(2α+1),w2=3ħ(ħ+1)tαsin⁡x2Γ(α+1)v3=-9ħ2(ħ+1)t2αsin⁡x2Γ(2α+1),w3=3ħ(ħ+1)2tαsin⁡x2Γ(α+1)+9ħ3t3αsin⁡x((-5+2cos⁡x)(Γ(α+1))2+2Γ(2α+1)sin2x)8(Γ(α+1))2Γ(3α+1),ψ(x,t)=∑n=07(vn+iwn).
The exact solution of (4.13) for α=1 is u=sinxe-3it/2. When ħ=-1, α=1,  the solution obtained by [13–16] is recovered as a special case.Example 4.3.
Consider the two dimensional NLFSE with trapping potential(4.15)iDαψ+12(∂2ψ∂x2+∂2ψ∂y2)+Vd(x,y)ψ+βd|ψ|2ψ=0,t≥0,0<α≤1,(x,y)∈[0,2π]×[0,2π].
Here Vd(x,y)=1-sin2xsin2y, βd=1, subject to the initial condition ψ(x,y,0)=sinxsiny. Solving the above equation we obtain the approximate solution in a series form
(4.16)v0=sin⁡xsin⁡y,w0=0v1=0,w1=2ħtαsin⁡xsin⁡yΓ(α+1)v2=4ħ2t2αsin⁡xsin⁡yΓ(2α+1),w2=2ħ(ħ+1)tαsin⁡xsin⁡yΓ(α+1)v3=-8ħ2(ħ+1)t2αsin⁡xsin⁡yΓ(2α+1)w3=-4ħ3t3αsin⁡xsin⁡y(2(1+sin2xsin2y)(Γ(α+1))2-Γ(2α+1)sin2xsin2y)(Γ(α+1))2Γ(3α+1)+2ħ(ħ+1)2tαsin⁡xsin⁡yΓ(α+1),ψ(x,y,t)=∑n=07(vn+iwn).
The exact solution of (4.15) for ħ=-1, α=1 is ψ=sinxsinye-2it.Example 4.4.
Consider the three-dimensional NLFSE with trapping potential(4.17)iDαψ+12(∂2ψ∂x2+∂2ψ∂y2+∂2ψ∂z2)-Vd(x,y,z)ψ-|ψ|2ψ=0,t≥0,0<α≤1,(4.18)(x,y,z)∈[0,2π]×[0,2π]×[0,2π],
where Vd(x,y,z)=1-sin2xsin2ysin2z subject to the initial condition
(4.19)ψ(x,y,z,0)=sin⁡xsin⁡ysin⁡z.
Solving the above equations,
(4.20)v0=sin⁡xsin⁡ysin⁡z,w0=0v1=0,w1=5ħtαsin⁡xsin⁡ysin⁡z2Γ(α+1)v2=-25ħ2t2αsin⁡xsin⁡ysin⁡z4Γ(2α+1),w2=-5ħ(ħ+1)tαsin⁡xsin⁡ysin⁡z2Γ(α+1).
Finally the approximate solution in the series form is
(4.21)ψ(x,y,z,t)=∑n=07(vn+iwn).
The exact solution of (4.17) for α=1 is ψ=sinxsinysinze-5it/2.


## 4.1. Nonlinear-Time Fractional Schrödinger Equation (NLFSE)

Settingψ(X,t)=v(X,t)+iw(X,t) and ψ(X,0)=v(X,0)+iw(X,0) in (2.1) leads to the following coupled system of equations(4.1)Dαv+12∇2w-Vdw-βd(wv2+w3)=0,(4.2)Dαw-12∇2v-Vdv+βd(v3+w2v)=0,subject to the initial conditions(4.3)v(X,0)=v0=F(X),w(X,0)=w0=G(X).Equations (4.1) and (4.3) suggest that we define the nonlinear operator as(4.4)NF1[φ1(X,t;q),φ2(X,t;q)]=Dαφ1(X,t;q)+12∇2φ2(X,t;q)-Vdφ2(X,t;q)-φ2(X,t;q)φ12(X,t;q)-φ23(X,t;q),(4.5)NF2[φ1(X,t;q),φ2(X,t;q)]=Dαφ2(X,t;q)-12∇2φ1(X,t;q)+Vdφ1(X,t;q)+φ13(X,t;q)+φ22(X,t;q)φ1(X,t;q),and the linear operator(4.6)£F[ϕ(X,t;q)]=Dα[ϕ(X,t;q)]with the property£F[c1(X)]=0£F[vm(X,t)-χmvm-1(X,t)]=ħR1,m[v⃗m-1,w⃗m-1](4.7)£F[wm(X,t)-χmwm-1(X,t)]=ħR2,m[v⃗m-1,w⃗m-1],where(4.8)R1,m[v⃗m-1,w⃗m-1]=Dαvm-1+12∇2wm-1-Vdwm-1-βd(∑k=0m-1vm-k-1∑l=0kvk-lwl+∑k=0m-1wm-k-1∑l=0kwk-lwl)R2,m[v⃗m-1,w⃗m-1]=Dαwm-1-12∇2vm-1+Vdvm-1+βd(∑k=0m-1wm-k-1∑l=0kwk-lvl+∑k=0m-1vm-k-1∑l=0kvk-lvl).Obviously, the solution of themth-order deformation (4.8) for m ≥ 1 becomes(4.9)vm(X,t)=χmvm-1(X,t)+ħJαR1,m[v⃗m-1,w⃗m-1]wm(X,t)=χmwm-1(X,t)+ħJαR2,m[v⃗m-1,w⃗m-1].Example 4.1.
We first consider the one-dimensional NLFSE with zero trapping potential (i.e.,Vd(x)=0) and βd=-1(4.10)iDtαψ(x,t)+12∂2ψ(x,t)∂x2+|ψ|2ψ=0,t≥0,0<α≤1.
subject to the initial condition ψ(x,0)=eix.
Solving the above equations(4.11)v0=cos⁡x,w0=sin⁡xv1=ħtαsin⁡x2Γ(α+1),w2=-ħtαcos⁡x2Γ(α+1)v2=ħ(ħ+1)tαsin⁡x2Γ(α+1)-ħ2t2αcos⁡x4Γ(2α+1),w2=-ħ(ħ+1)tαcos⁡x2Γ(α+1)-ħ2t2αsin⁡x4Γ(2α+1)v3=(ħ(ħ+1)2tα2Γ(α+1)-ħ3t3α(5(Γ(α+1))2-2Γ(2α+1))8(Γ(α+1))2Γ(3α+1))sin⁡x-ħ2(ħ+1)t2αcos⁡x2Γ(2α+1),w3=(-ħ(ħ+1)2tα2Γ(α+1)+ħ3t3α(5(Γ(α+1))2-2Γ(2α+1))8(Γ(α+1))2Γ(3α+1))cos⁡x-ħ2(ħ+1)t2αsin⁡x2Γ(2α+1),and so on, in this manner the rest of the components can be obtained. Therefore, the approximate solution is
(4.12)ψ(x,t)=∑n=07(vn+iwn).
The exact solution of (4.10) for α=1 is ψ=ei(x+(t/2)). When ħ=-1, α=1, the solution obtained by [13–16] is recovered as a special case.Example 4.2.
Consider the one-dimensional NLFSE with trapping potential, that is,Vd(x)=cos2x and βd=-1(4.13)iDtαψ(x,t)+12∂2ψ(x,t)∂x2-ψcos⁡2⁡x-|ψ|2ψ=0,t≥0,0<α≤1.
subject to the initial condition ψ(x,0)=sinx.
Solving above equations we obtain the solution in a series form(4.14)v0=sin⁡x,w0=0v1=0,w1=3ħtαsin⁡x2Γ(α+1)v2=-9ħ2t2αsin⁡x4Γ(2α+1),w2=3ħ(ħ+1)tαsin⁡x2Γ(α+1)v3=-9ħ2(ħ+1)t2αsin⁡x2Γ(2α+1),w3=3ħ(ħ+1)2tαsin⁡x2Γ(α+1)+9ħ3t3αsin⁡x((-5+2cos⁡x)(Γ(α+1))2+2Γ(2α+1)sin2x)8(Γ(α+1))2Γ(3α+1),ψ(x,t)=∑n=07(vn+iwn).
The exact solution of (4.13) for α=1 is u=sinxe-3it/2. When ħ=-1, α=1,  the solution obtained by [13–16] is recovered as a special case.Example 4.3.
Consider the two dimensional NLFSE with trapping potential(4.15)iDαψ+12(∂2ψ∂x2+∂2ψ∂y2)+Vd(x,y)ψ+βd|ψ|2ψ=0,t≥0,0<α≤1,(x,y)∈[0,2π]×[0,2π].
Here Vd(x,y)=1-sin2xsin2y, βd=1, subject to the initial condition ψ(x,y,0)=sinxsiny. Solving the above equation we obtain the approximate solution in a series form
(4.16)v0=sin⁡xsin⁡y,w0=0v1=0,w1=2ħtαsin⁡xsin⁡yΓ(α+1)v2=4ħ2t2αsin⁡xsin⁡yΓ(2α+1),w2=2ħ(ħ+1)tαsin⁡xsin⁡yΓ(α+1)v3=-8ħ2(ħ+1)t2αsin⁡xsin⁡yΓ(2α+1)w3=-4ħ3t3αsin⁡xsin⁡y(2(1+sin2xsin2y)(Γ(α+1))2-Γ(2α+1)sin2xsin2y)(Γ(α+1))2Γ(3α+1)+2ħ(ħ+1)2tαsin⁡xsin⁡yΓ(α+1),ψ(x,y,t)=∑n=07(vn+iwn).
The exact solution of (4.15) for ħ=-1, α=1 is ψ=sinxsinye-2it.Example 4.4.
Consider the three-dimensional NLFSE with trapping potential(4.17)iDαψ+12(∂2ψ∂x2+∂2ψ∂y2+∂2ψ∂z2)-Vd(x,y,z)ψ-|ψ|2ψ=0,t≥0,0<α≤1,(4.18)(x,y,z)∈[0,2π]×[0,2π]×[0,2π],
where Vd(x,y,z)=1-sin2xsin2ysin2z subject to the initial condition
(4.19)ψ(x,y,z,0)=sin⁡xsin⁡ysin⁡z.
Solving the above equations,
(4.20)v0=sin⁡xsin⁡ysin⁡z,w0=0v1=0,w1=5ħtαsin⁡xsin⁡ysin⁡z2Γ(α+1)v2=-25ħ2t2αsin⁡xsin⁡ysin⁡z4Γ(2α+1),w2=-5ħ(ħ+1)tαsin⁡xsin⁡ysin⁡z2Γ(α+1).
Finally the approximate solution in the series form is
(4.21)ψ(x,y,z,t)=∑n=07(vn+iwn).
The exact solution of (4.17) for α=1 is ψ=sinxsinysinze-5it/2.


## 5. Closing Comments

The basic goal of this work has been to extend the works made in the nonlinear physical problem to construct solutions for nonlinear Schrödinger equation with time-fractional derivatives. The goal has been achieved and new solutions of have been derived for nonlinear equations with time fractional derivatives. The proposed approach works successfully in handling nonlinear fractional Schrödinger equations directly with a minimum size of calculations.Borhanifer and Abazari [11] applied the differential transformation method (DTM) to solving Schrödinger and coupled Schrödinger equations. The major lacks of DTM are that it requires transformation, and the given differential equation and related initial conditions are transformed into a recurrence equation that finally leads to the solution of a system of algebraic equations as coefficients of a power series solution. The main disadvantage of DTM is that it requires transformation, which will be complicated and computational cost will be too much. This emphasizes the fact that the presented approach can be used in a wider class of system of nonlinear fractional differential equations. HAM is a powerful and efficient technique in finding exact and approximate solutions for linear and nonlinear models. HAM provides more realistic solutions that converge very rapidly in real physical problems. The numerical examples show that the solutions are in good agreement with their respective exact solutions for α=1.In the last, we present theħ-curves to see the convergent region of ħ as in Figures 1–4. We plot the imaginary and real part of each example in the same figure, and the convergent region is the region of intersection between the convergent regions of the imaginary and real parts. In all figures ħ=-1 is in the convergent region.Figure 1
Theħ-curve for Example 4.1  x=0.2, t=0.2, α=0.7.Figure 2
Theħ-curve for Example 4.2  x=0.2, t=0.2,α=0.7.Figure 3
Theħ-curve for Example 4.3  x=0.1, y=0.1, t=0.5, α=0.9.Figure 4
Theħ-curve for Example 4.4  x=0.1, y=0.1, z=0.1, t=0.5, α=0.9.



---

*Source: 197068-2012-03-28.xml*