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.

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 [15] 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) [611] 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) [1720].

The aim of this paper is to solve the Schrödinger equation with fractional order using the homotopy analysis method. By introducing an embedding parameter 𝑞 the nonlinear fractional differential equation is converted to a linear fractional differential equation at 𝑞=0. When 𝑞 evolves, the differential equation becomes the original one at 𝑞=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: 𝑖𝜕𝛼𝜓(𝑋,𝑡)𝜕𝑡𝛼+122𝜓(𝑋,𝑡)𝑉𝑑(𝑋)𝜓𝛽𝑑||𝜓||2𝜓=0;𝑡0,𝜓(𝑋,0)=𝜓0(𝑋),𝑋𝑑,(2.1) where 𝑉𝑑(𝑋) is the trapping potential and 𝛽𝑑 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 𝐿1[𝑎,𝑏] is given by 𝐽𝜇1𝑓(𝑥)=Γ(𝜇)𝑥0(𝑥𝑡)𝜇1𝐽𝑓(𝑡)𝑑𝑡;𝜇>0,(2.2)0𝑓(𝑥)=𝑓(𝑥).(2.3) This integral operator has the following properties:(i)𝐽𝛼𝐽𝛽=𝐽𝛼+𝛽=𝐽𝛽𝐽𝛼,𝛼,𝛽>0,(ii)𝐽𝛼(𝑥𝑎)𝛾=Γ(𝛾+1)Γ(𝛼+𝛾+1)(𝑥𝑎)𝛼+𝛾,𝛼>0,𝛾>1.

Definition 2.2. The Caputo definition of fractal derivative operator is given by 𝐷𝜇1𝑓(𝑥)=Γ(𝑚𝜇)𝑡0(𝑥𝜏)𝑚𝜇1𝑓(𝑚)(𝜏)𝑑𝜏𝑚1<𝜇𝑚,𝑚𝑁,𝑥>0.(2.4) It has the following two basic properties: 𝐷𝜇𝐽𝜇𝑓(𝑥)=𝑓(𝑥),(2.5) and 𝐽𝜇𝐷𝜇𝑓(𝑥)=𝑓(𝑥)𝑚1𝑘=0𝑓(𝑘)0+(𝑥)𝑘𝑘!,𝑥>0.(2.6) 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:𝒩[]𝜓(𝑋,𝑡;𝑞)=0,(3.1)

where 𝒩 is a nonlinear fractional operator, 𝑋 and 𝑡 denote the independent variables, and 𝜓 is an unknown function. By means of the HAM, we first construct the so-called zeroth-order deformation equation(1𝑞)£𝜙(𝑋,𝑡;𝑞)𝜓0(𝑋,𝑡)=𝑞𝒩[𝜙](𝑋,𝑡;𝑞),(3.2)

where 𝑞[0,1] is the embedding parameter, 0 is an auxiliary parameter, £ is an auxiliary linear operator, 𝜙(𝑋,𝑡;𝑞) are unknown functions, and 𝜓0(𝑋,𝑡) are initial guesses of 𝜙(𝑋,𝑡;𝑞). It is obvious that for 𝑞=0 and 𝑞=1, (3.3) becomes𝜙(𝑋,𝑡;0)=𝜓0(𝑋,𝑡),𝜙(𝑋,𝑡;1)=𝜓(𝑋,𝑡),(3.3)

respectively. Thus as 𝑞=0 increases from 0 to 1, the solution 𝜙(𝑋,𝑡;𝑞) varies from the initial guess 𝜓0(𝑋,𝑡) to the exact solution 𝜓(𝑋,𝑡). Expanding 𝜙(𝑋,𝑡;𝑞) in Taylor series with respect to 𝑞, one has𝜙(𝑋,𝑡;𝑞)=𝜓0(𝑋,𝑡)++𝑚=1𝜓𝑚(𝑋,𝑡)𝑞𝑚,(3.4)

where𝜓𝑚1(𝑋,𝑡)=𝜕𝑚!𝑚𝜙(𝑋,𝑡;𝑞)𝜕𝑞𝑚||||𝑞=0.(3.5)

The convergence of the series (3.4) depends upon the auxiliary parameter . If it is convergent at 𝑞=1, one has𝜓(𝑋,𝑡)=𝜓0(𝑋,𝑡)++𝑚=1𝜓𝑚(𝑋,𝑡),(3.6)

which must be one of the solutions of the original nonlinear equations, as proved by Liao [17]. Define the vectors𝜓𝑛=𝜓0(𝑋,𝑡),𝜓1(𝑋,𝑡),𝜓2(𝑋,𝑡),,𝜓𝑛.(𝑋,𝑡)(3.7) Differentiating the zeroth-order deformation (3.2) 𝑚-times with respect to 𝑞 and then dividing them by 𝑚!, and finally setting 𝑞=0, we get the following 𝑚th-order deformation equation:£𝜓𝑚(𝑋,𝑡)𝜒𝑚𝜓𝑚1(𝑋,𝑡)=𝑅𝑚𝜓𝑚1,(3.8)

where𝑅𝑚𝜓𝑚1=1𝜕(𝑚1)!𝑚1𝒩[𝜙](𝑋,𝑡;𝑞)𝜕𝑞𝑚1||||𝑞=0,𝜒𝑚=0;𝑚=1,1;𝑚>1.(3.9)   

It should be emphasized that 𝜓𝑚(𝑋,𝑡;𝑞) for 𝑚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 𝜓(𝑋,𝑡)=𝑣(𝑋,𝑡)+𝑖w(𝑋,𝑡) and 𝜓(𝑋,0)=𝑣(𝑋,0)+𝑖𝑤(𝑋,0) in (2.1) leads to the following coupled system of equations𝐷𝛼1𝑣+22𝑤𝑉𝑑𝑤𝛽𝑑𝑤𝑣2+𝑤3𝐷=0,(4.1)𝛼1𝑤22𝑣𝑉𝑑𝑣+𝛽𝑑𝑣3+𝑤2𝑣=0,(4.2)

subject to the initial conditions𝑣(𝑋,0)=𝑣0=𝐹(𝑋),𝑤(𝑋,0)=𝑤0=𝐺(𝑋).(4.3)

Equations (4.1) and (4.3) suggest that we define the nonlinear operator as𝒩1𝜑1(𝑋,𝑡;𝑞),𝜑2(𝑋,𝑡;𝑞)=𝐷𝛼𝜑11(𝑋,𝑡;𝑞)+22𝜑2(𝑋,𝑡;𝑞)𝑉𝑑𝜑2(𝑋,𝑡;𝑞)𝜑2(𝑋,𝑡;𝑞)𝜑21(𝑋,𝑡;𝑞)𝜑32𝒩(𝑋,𝑡;𝑞),(4.4)2𝜑1(𝑋,𝑡;𝑞),𝜑2(𝑋,𝑡;𝑞)=𝐷𝛼𝜑21(𝑋,𝑡;𝑞)22𝜑1(𝑋,𝑡;𝑞)+𝑉𝑑𝜑1(𝑋,𝑡;𝑞)+𝜑31(𝑋,𝑡;𝑞)+𝜑22(𝑋,𝑡;𝑞)𝜑1(𝑋,𝑡;𝑞),(4.5)

and the linear operator£[𝜙](𝑋,𝑡;𝑞)=𝐷𝛼[𝜙](𝑋,𝑡;𝑞)withtheproperty£𝑐1£(𝑋)=0𝑣𝑚(𝑋,𝑡)𝜒𝑚𝑣𝑚1(𝑋,𝑡)=𝑅1,𝑚𝑣𝑚1,𝑤𝑚1£(4.6)𝑤𝑚(𝑋,𝑡)𝜒𝑚𝑤𝑚1(𝑋,𝑡)=𝑅2,𝑚𝑣𝑚1,𝑤𝑚1,(4.7)

where𝑅1,𝑚𝑣𝑚1,𝑤𝑚1=𝐷𝛼𝑣𝑚1+122𝑤𝑚1𝑉𝑑𝑤𝑚1𝛽𝑑𝑚1𝑘=0𝑣𝑘𝑚𝑘1𝑙=0𝑣𝑘𝑙𝑤𝑙+𝑚1𝑘=0𝑤𝑘𝑚𝑘1𝑙=0𝑤𝑘𝑙𝑤𝑙𝑅2,𝑚𝑣𝑚1,𝑤𝑚1=𝐷𝛼𝑤𝑚1122𝑣𝑚1+𝑉𝑑𝑣𝑚1+𝛽𝑑𝑚1𝑘=0𝑤𝑘𝑚𝑘1𝑙=0𝑤𝑘𝑙𝑣𝑙+𝑚1𝑘=0𝑣𝑘𝑚𝑘1𝑙=0𝑣𝑘𝑙𝑣𝑙.(4.8)

Obviously, the solution of the 𝑚th-order deformation (4.8) for m ≥ 1 becomes𝑣𝑚(𝑋,𝑡)=𝜒𝑚𝑣𝑚1(𝑋,𝑡)+𝐽𝛼𝑅1,𝑚𝑣𝑚1,𝑤𝑚1𝑤𝑚(𝑋,𝑡)=𝜒𝑚𝑤𝑚1(𝑋,𝑡)+𝐽𝛼𝑅2,𝑚𝑣𝑚1,𝑤𝑚1.(4.9)

Example 4.1. We first consider the one-dimensional NLFSE with zero trapping potential (i.e.,𝑉𝑑(𝑥)=0) and 𝛽𝑑=1𝑖𝐷𝛼𝑡1𝜓(𝑥,𝑡)+2𝜕2𝜓(𝑥,𝑡)𝜕𝑥2+||𝜓||2𝜓=0,𝑡0,0<𝛼1.(4.10) subject to the initial condition 𝜓(𝑥,0)=𝑒𝑖𝑥.
Solving the above equations 𝑣0=cos𝑥,𝑤0𝑣=sin𝑥1=𝑡𝛼sin𝑥2Γ(𝛼+1),𝑤2=𝑡𝛼cos𝑥𝑣2Γ(𝛼+1)2=(+1)𝑡𝛼sin𝑥2Γ(𝛼+1)2𝑡2𝛼cos𝑥4Γ(2𝛼+1),𝑤2=(+1)𝑡𝛼cos𝑥2Γ(𝛼+1)2𝑡2𝛼sin𝑥𝑣4Γ(2𝛼+1)3=(+1)2𝑡𝛼2Γ(𝛼+1)3𝑡3𝛼5(Γ(𝛼+1))22Γ(2𝛼+1)8(Γ(𝛼+1))2Γ(3𝛼+1)sin𝑥2(+1)𝑡2𝛼cos𝑥,𝑤2Γ(2𝛼+1)3=(+1)2𝑡𝛼+2Γ(𝛼+1)3𝑡3𝛼5(Γ(𝛼+1))22Γ(2𝛼+1)8(Γ(𝛼+1))2Γ(3𝛼+1)cos𝑥2(+1)𝑡2𝛼sin𝑥,2Γ(2𝛼+1)(4.11)and so on, in this manner the rest of the components can be obtained. Therefore, the approximate solution is 𝜓(𝑥,𝑡)=7𝑛=0𝑣𝑛+𝑖𝑤𝑛.(4.12) The exact solution of (4.10) for 𝜶=1 is 𝝍=𝑒𝑖(𝑥+(𝑡/2)). When =1, 𝜶=1, the solution obtained by [1316] is recovered as a special case.

Example 4.2. Consider the one-dimensional NLFSE with trapping potential, that is, 𝑉𝑑(𝑥)=cos2𝑥 and 𝛽𝑑=1𝑖𝐷𝛼𝑡1𝜓(𝑥,𝑡)+2𝜕2𝜓(𝑥,𝑡)𝜕𝑥2𝜓cos2||𝜓||𝑥2𝜓=0,𝑡0,0<𝛼1.(4.13) subject to the initial condition 𝜓(𝑥,0)=sin𝑥.
Solving above equations we obtain the solution in a series form 𝑣0=sin𝑥,𝑤0𝑣=01=0,𝑤1=3𝑡𝛼sin𝑥𝑣2Γ(𝛼+1)2=92𝑡2𝛼sin𝑥4Γ(2𝛼+1),𝑤2=3(+1)𝑡𝛼sin𝑥𝑣2Γ(𝛼+1)3=92(+1)𝑡2𝛼sin𝑥,𝑤2Γ(2𝛼+1)3=3(+1)2𝑡𝛼sin𝑥+2Γ(𝛼+1)93𝑡3𝛼sin𝑥(5+2cos𝑥)(Γ(𝛼+1))2+2Γ(2𝛼+1)sin2𝑥8(Γ(𝛼+1))2Γ,(3𝛼+1)𝜓(𝑥,𝑡)=7𝑛=0𝑣𝑛+𝑖𝑤𝑛.(4.14) The exact solution of (4.13) for 𝛼=1 is 𝑢=𝑠in𝑥𝑒3𝑖𝑡/2. When =1, 𝛼=1, the solution obtained by [1316] is recovered as a special case.

Example 4.3. Consider the two dimensional NLFSE with trapping potential 𝑖𝐷𝛼1𝜓+2𝜕2𝜓𝜕𝑥2+𝜕2𝜓𝜕𝑦2+𝑉𝑑(𝑥,𝑦)𝜓+𝛽𝑑||𝜓||2[]×[].𝜓=0,𝑡0,0<𝛼1,(𝑥,𝑦)0,2𝜋0,2𝜋(4.15) Here 𝑉𝑑(𝑥,𝑦)=1sin2𝑥sin2𝑦, 𝛽𝑑=1, subject to the initial condition 𝜓(𝑥,𝑦,0)=sin𝑥sin𝑦. Solving the above equation we obtain the approximate solution in a series form 𝑣0=sin𝑥sin𝑦,𝑤0𝑣=01=0,𝑤1=2𝑡𝛼sin𝑥sin𝑦𝑣Γ(𝛼+1)2=42𝑡2𝛼sin𝑥sin𝑦Γ(2𝛼+1),𝑤2=2(+1)𝑡𝛼sin𝑥sin𝑦Γ𝑣(𝛼+1)3=82(+1)𝑡2𝛼sin𝑥sin𝑦𝑤Γ(2𝛼+1)3=43𝑡3𝛼2sin𝑥sin𝑦1+sin2𝑥sin2𝑦(Γ(𝛼+1))2Γ(2𝛼+1)sin2𝑥sin2𝑦(Γ(𝛼+1))2+Γ(3𝛼+1)2(+1)2𝑡𝛼sin𝑥sin𝑦,𝜓Γ(𝛼+1)(𝑥,𝑦,𝑡)=7𝑛=0𝑣𝑛+𝑖𝑤𝑛.(4.16) The exact solution of (4.15) for =1, 𝛼=1 is 𝜓=sin𝑥sin𝑦𝑒2𝑖𝑡.

Example 4.4. Consider the three-dimensional NLFSE with trapping potential 𝑖𝐷𝛼1𝜓+2𝜕2𝜓𝜕𝑥2+𝜕2𝜓𝜕𝑦2+𝜕2𝜓𝜕𝑧2𝑉𝑑||𝜓||(𝑥,𝑦,𝑧)𝜓2[]×[]×[]𝜓=0,𝑡0,0<𝛼1,(4.17)(𝑥,𝑦,𝑧)0,2𝜋0,2𝜋0,2𝜋,(4.18) where 𝑉𝑑(𝑥,𝑦,𝑧)=1sin2𝑥sin2𝑦sin2𝑧 subject to the initial condition 𝜓(𝑥,𝑦,𝑧,0)=sin𝑥sin𝑦sin𝑧.(4.19) Solving the above equations, 𝑣0=sin𝑥sin𝑦sin𝑧,𝑤0𝑣=01=0,𝑤1=5𝑡𝛼sin𝑥sin𝑦sin𝑧𝑣2Γ(𝛼+1)2=252𝑡2𝛼sin𝑥sin𝑦sin𝑧4Γ(2𝛼+1),𝑤2=5(+1)𝑡𝛼sin𝑥sin𝑦sin𝑧.2Γ(𝛼+1)(4.20) Finally the approximate solution in the series form is 𝜓(𝑥,𝑦,𝑧,𝑡)=7𝑛=0𝑣𝑛+𝑖𝑤𝑛.(4.21) The exact solution of (4.17) for 𝛼=1 is 𝜓=sin𝑥sin𝑦sin𝑧𝑒5𝑖𝑡/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 14. 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.

Acknowledgments

The author N. A. Khan is thankful and grateful to the dean of the Faculty of Sciences, University of Karachi, Pakistan for supporting and facilitating this research work. The author M. Jamil is highly thankful and grateful to the Abdul Salam School of Mathematical Sciences, GC University, Lahore, Pakistan, Department of Mathematics, NED University of Engineering and Technology, Pakistan, and also Higher Education Commission of Pakistan for generous support and facilitating this research work.