Time Fractional Schrodinger Equation Revisited

Abstract

The time fractional Schrodinger equation (TFSE) for a nonrelativistic particle is derived on the basis of the Feynman path integral method by extending it initially to the case of a “free particle” obeying fractional dynamics, obtained by replacing the integer order derivatives with respect to time by those of fractional order. The equations of motion contain quantities which have “fractional” dimensions, chosen such that the “energy” has the correct dimension [𝑀𝐿²/𝑇²]. The action 𝑆 is defined as a fractional time integral of the Lagrangian, and a “fractional Planck constant” is introduced. The TFSE corresponds to a “subdiffusion” equation with an imaginary fractional diffusion constant and reproduces the regular Schrodinger equation in the limit of integer order. The present work corrects a number of errors in Naber’s work. The correct continuity equation for the probability density is derived and a Green function solution for the case of a “free particle” is obtained. The total probability for a “free” particle is shown to go to zero in the limit of infinite time, in contrast with Naber’s result of a total probability greater than unity. A generalization to the case of a particle moving in a potential is also given.

1. Introduction

There has been an explosive research output in recent years in the application of methods of fractional calculus [1–13] to the study of quantum phenomena [14–42]. The well-known Schrodinger equation with a first-order derivative in time and second-order derivatives in space coordinates was given by Schrodinger as an Ansatz. The Schrodinger equation has been generalized to (i) a space fractional Schrodinger equation involving noninteger order space derivatives but retaining first-order time derivative [14–18], (ii) a time fractional Schrodinger equation involving non-integer order time derivative but retaining the second-order space derivatives [19], or (iii) more general fractional Schrodinger equation where both time and space derivatives are of non-integer order [20–26]. The fractional Schrodinger equation has also been obtained by using a fractional generalization of the Laplacian operator [20] and by using a fractional variation principle and a fractional Klein-Gordon equation [36]. In all these cases the fractional derivatives employed have been the regular fractional derivatives of the Riemann-Liouville type or the Caputo type (generally used in physical applications with initial conditions) which are both nonlocal in nature. The fractional derivative which is nonlocal by definition can be made “local” by a limiting process as shown by Kolwankar and Gangal [41]. Highly irregular and nowhere differentiable functions can be analyzed locally using these local fractional derivatives. The Heisenberg principle in the fractional context has been investigated using local fractional Fourier analysis [42].

The Schrodinger equation for a free particle has the appearance of a diffusion equation with an imaginary diffusion coefficient. This suggests a method of deriving the Schrodinger equation as has been done using the Feynman path integral technique [43–45] based on the Gaussian probability distribution in the space of all possible paths. In other words, the classical Brownian motion leads to the Schrodinger equation in quantum mechanics. As far as deriving the fractional Schrodinger equation is concerned, the path integral approach for the Brownian-like paths for the Levy stable processes which leads to the classical space fractional diffusion equation has been extended to the Levy-like quantum paths leading to the space fractional Schrodinger equation (SFSE) in the seminal papers of Laskin [14–18]. It may be noted that in this case, the time derivative is still the integer first-order derivative; only the space part is of fractional order. The SFSE still retains the Markovian character and other fundamental aspects such as the Hermiticity of the Hamilton operator. Parity conservation and the current density have been explored in the space fractional quantum mechanics in terms of the Riesz fractional derivative. Applications of SFSE cover the dynamics of a free particle, particle in an infinite potential well, fractional Bohr atom, and the quantum fractional oscillator. Thus the space fractional Schrodinger equation appears to have been well established [18]. This theory has been further generalized recently within the frame work of tempered ultradistributions [39]. Thus the theory of SFSE can be considered fully established from the point of view of the Feynman path integral technique.

The fractional time derivative was introduced into the Schrodinger equation by Naber [19] by simply replacing the first-order time derivative by a derivative of non-integer order and retaining the second-order space derivatives intact. The resulting equation is referred to as the time fractional Schrodinger equation (TFSE). He did not derive the TFSE using the path integral or any other method. Naber carried out the time fractional modification to the Schrodinger equation in analogy with time fractional diffusion equation [19] but included the imaginary number 𝑖 raised to a fractional power (the fractional degree being the same as the fractional order of the time derivative), implying a sort of the Wick rotation. In Naber’s opinion [19], the TFSE is equivalent to the usual Schrodinger equation, but with a time-dependent Hamiltonian. He obtained the solutions for a free particle and a particle in a potential well. A lot of subsequent work has been done on the TFSE, mostly based on Naber’s work [21–23, 28], including its generalization into space-time fractional quantum dynamics by including non-integer order derivatives in both time and space. Yet some basic questions have not been addressed. It has been observed that TFSE describes non-Markovian evolution and that the Green function in the form of the Mittag-Leffler function does not satisfy Stone’s theorem on one-parameter unitary groups and the semiclassical approximation in terms of the classical action is not defined [31]. There has been no derivation of TFSE on a basis similar to that of SFSE and it is the purpose of the present paper to rectify this lacuna. Since the path integral method of deriving the space fractional part of the Schrodinger equation is well established, the present paper concentrates only on deriving the time fractional Schrodinger equation from the Feynman path integral approach, leading to the time fractional Schrodinger equation as given by Naber. It may be pointed out that some results of Naber, such as the total probability being greater than unity, are difficult to understand physically. Moreover, several major errors in Naber’s paper have gone unnoticed and in fact have been repeated by workers who have followed his work. Furthermore, some of these authors have introduced errors of their own. Since many of the conclusions in Naber’s paper are based on derivations which include these errors, it calls for a reexamination of Naber’s generalization to the time fractional Schrodinger equation.

The present paper derives the time fractional Schrodinger equation using the Feynman path integral technique. It concentrates on the time fractional part only and not on the space fractional part as the theory of the latter has been well established in the works of Laskin [14–18]. Furthermore, the paper considers only the Caputo-type nonlocal fractional derivatives and not the local fractional derivatives discussed earlier. The paper starts from a generalization of the classical dynamics into fractional dynamics of a free particle and then adapting the Feynman technique derives the correct equations for TFSE. It is demonstrated that Naber’s result of probability being greater than unity is spurious and is a result of the ad hoc raising of the imaginary number 𝑖 to a fractional power. The correct continuity equation for the probability density is also derived. The paper concludes with some new results.

2. Feynman Path Integral Method

The starting point for the Feynman method [43–46] is the classical Lagrangian 𝐿=𝐿(𝑥,𝑥̇,𝑡) and the action 𝑆=∫𝐿(𝑥,𝑥̇,𝑡)𝑑𝑡 constructed from it. However, in view of the generalization to fractional calculus methods to be carried out later, the equations of motion of a classical particle in one dimension in the usual notation are considered first: 𝑚𝑑𝑥𝑑𝑡=𝑝,𝑑𝑝𝑑𝑡=𝐹. (1) Integrating with respect to time yields 𝑥=𝑥₀+1𝑚^𝑡₀𝑝(𝜏)𝑑𝜏,𝑝=𝑝₀+^𝑡₀𝐹(𝜏)𝑑𝜏. (2)(3) In the usual notation the Lagrangian is given by 𝐿(𝑥,𝑥̇,𝑡)=𝑇−𝑉 (4) and the action is given by 𝑆=^𝑡₀𝐿(𝑥,𝑥̇,𝜏)𝑑𝜏. (5) An outline of the Feynman path integral method is presented following very closely the account given by Feynman and Hibbs [43]. The essence of the Feynman path integral approach to quantum mechanics is in the probability amplitude (also known as the propagator or the Green function) 𝐾(𝑥_𝑏,𝑡_𝑏;𝑥_𝑎,𝑡_𝑎) for a particle starting from a position 𝑥_𝑎 at time 𝑡_𝑎 to reach a position 𝑥_𝑏 at a later time 𝑡_𝑏, which arises from the contributions from all trajectories from 𝑥_𝑎 to 𝑥_𝑏: 𝐾𝑥_𝑏,𝑡_𝑏;𝑥_𝑎,𝑡_𝑎=_allpaths𝜙[𝑥(𝑡)], (6) where the contribution from each of the paths has the form 𝜙[𝑥(𝑡)]=const.exp𝑖𝑆ℏ. (7) Here 𝑆 is the action defined in (5) and ℏ is Planck’s constant, the quantum of action. The time integral of the Lagrangian is to be taken along the path in question. Restricting to one dimension, the probability amplitude can be written as 𝐾𝑥_𝑏,𝑡_𝑏;𝑥_𝑎,𝑡_𝑎=^𝑏_𝑎𝑖expℏ^𝑡𝑏_𝑡𝑎𝐿𝑑𝑡𝔇𝑥(𝑡). (8) The symbol 𝔇 indicates the fact that the operation of integration is carried over all paths from 𝑎 to 𝑏.

The wave function 𝜓(𝑥_𝑏,𝑡_𝑏) gives the total probability amplitude to arrive at 𝑥_𝑏 at 𝑡_𝑏 satisfying (9), where the integral is taken over all possible values of 𝑥_𝑎𝜓𝑥_𝑏,𝑡_𝑏=^∞_−∞𝐾𝑥_𝑏,𝑡_𝑏;𝑥_𝑎,𝑡_𝑎𝜓𝑥_𝑎,𝑡_𝑎𝑑𝑥_𝑎. (9) The kernel 𝐾 can be computed by first carrying out a “time slicing” operation by dividing the time interval from 𝑡_𝑎 to 𝑡_𝑏 into 𝑁 segments of duration 𝑡𝜀=^𝑏−𝑡_𝑎𝑁, (10) where 𝑡_𝑎=𝑡₀<𝑡₁<𝑡₂<⋯<𝑡_𝑁−1<𝑡_𝑁=𝑡_𝑏;𝐾𝑥_𝑏,𝑡_𝑏;𝑥_𝑎,𝑡_𝑎=lim_𝜀→01𝐴⋯[]exp𝑖𝑆𝑏,𝑎ℏ×𝑑𝑥₁𝐴𝑑𝑥₂𝐴𝑑𝑥₃𝐴⋯𝑑𝑥_𝑁−1𝐴, (11)(12) where 𝑆[𝑏,𝑎]=^𝑡𝑏_𝑡𝑎𝐿(𝑥,𝑥̇,𝑡)𝑑𝑡 (13) is a line integral taken over the trajectory passing through the point 𝑥(𝑡). The constant 𝐴 is a normalizing factor.

The Schrodinger equation for a free particle in one dimension is derived by considering a special case of (9), which describes the evolution of the wave function from a time 𝑡_𝑎 to a time 𝑡_𝑏, when 𝑡_𝑏 differs from 𝑡_𝑎 by an infinitesimal amount 𝜀 and applied to the case of a free particle. This step is based on the fact that the semiclassical approximation is valid not only in the limit of ℏ→0 but also in the limit of small time interval [45]. The Kernel is proportional to the exponential of (𝑖/ℏ) times the classical action for the infinitesimal time interval 𝜀=𝑡_𝑏−𝑡_𝑎. With an obvious change of notation 𝑥_𝑏=𝑥, 𝑥_𝑎=𝑥₀, 𝑡_𝑎=𝑡, 𝑡_𝑏=𝑡+𝜀 and using the fact that the particle is free, (2) yields 𝑝₀=𝑚𝑥−𝑥₀𝜀 (14) and (4) yields 𝑚𝐿=𝑇=𝑥−𝑥₀²2𝜀² (15) and (5) yields for the action 𝑚𝑆=𝜀𝐿=𝑥−𝑥₀².2𝜀 (16) As a consequence, (9) becomes 𝜓(𝑥,𝑡+𝜀)=^∞_−∞1𝐴exp𝑖𝑚𝑥−𝑥₀²𝜓𝑥2𝜀⁰,𝑡𝑑𝑥₀. (17) If 𝑥 differs appreciably from 𝑥₀, the exponential in (17) oscillates very rapidly and the integral over 𝑥₀ contributes a very small value and only those paths which are very close to 𝑥 give significant contributions. Changing the variable in the integral from 𝑥₀ to 𝜂=𝑥−𝑥₀ makes it 𝜓(𝑥₀,𝑡)=𝜓(𝑥+𝜂,𝑡). Since both 𝜀 and 𝜂 are small quantities, 𝜓(𝑥,𝑡+𝜀) may be expanded in Taylor’s series and only up to terms of order 𝜀 are retained. On the right-hand side, 𝜓(𝑥+𝜂,𝑡) may be expanded in Taylor’s series in powers of 𝜂, retaining terms up to second order in 𝜂 (the integral involving the first-order term vanishes). Then (17) becomes 𝜓(𝑥,𝑡)+𝜀𝜕𝜓=𝜕𝑡_∞−∞1𝐴exp−𝑚𝜂²×𝜓()2𝑖ℏ𝜀𝑥,𝑡+𝜂𝜕𝜓+𝜂𝜕𝑥2²𝜕²𝜓𝜕𝑥²𝑑𝜂. (18) On the right-hand side the middle term vanishes on integration. It follows by equating the leading terms on both sides 𝜓(𝑥,𝑡)=^∞_−∞1𝐴exp−𝑚𝜂²𝜓()2𝑖ℏ𝜀𝑥,𝑡𝑑𝜂 (19) Hence𝐴=^∞_−∞exp−𝑚𝜂²2𝑖ℏ𝜀𝑑𝜂=2𝜋𝑖ℏ𝜀𝑚, (20a)^∞_−∞1𝐴exp−𝑚𝜂²𝜂2𝑖ℏ𝜀2²𝜕²𝜓𝜕𝑥²𝑑𝜂=𝜀𝑖ℏ𝜕2𝑚²𝜓𝜕𝑥². (20b)Equating the remaining terms results in 𝜕𝜓=𝜕𝑡𝑖ℏ𝜕2𝑚²𝜓𝜕𝑥². (21) This can be recognized as the diffusion equation with an imaginary diffusion coefficient or the Schrodinger equation for a free particle in quantum mechanics.

These considerations can be easily extended to the case of a particle moving in a potential field by incorporating a potential term 𝑉(𝑥,𝑡) in the Lagrangian 𝐿=𝑇−𝑉 in (15). This will necessitate [43] incorporating an additional factor {1−(𝑖𝜀/ℏ)𝑉(𝑥,𝑡)} in (18), which becomes 𝜓(𝑥,𝑡)+𝜀𝜕𝜓=𝜕𝑡^∞_−∞1𝐴exp−𝑚𝜂²2𝑖ℏ𝜀1−𝑖𝜀ℏ𝑉(×𝜓()𝑥,𝑡)𝑥,𝑡+𝜂𝜕𝜓+𝜂𝜕𝑥2²𝜕²𝜓𝜕𝑥²𝑑𝜂. (22) Then (21) becomes 𝜕𝜓=𝜕𝑡𝑖ℏ𝜕2𝑚²𝜓𝜕𝑥²−𝑖ℏ𝑉𝜓. (23) Multiplying both sides by −ℏ/𝑖 results in the standard Schrodinger equation of quantum mechanics: 𝑖ℏ𝜕𝜓ℏ𝜕𝑡=−²𝜕2𝑚²𝜓𝜕𝑥²+𝑉𝜓. (24) It is to be noted that the imaginary number 𝑖 on the left-hand side is not arbitrary; it arises from the coefficient of the potential term 𝑉 in (23) but ultimately from the coefficient of 𝑆/ℏ in (7). This minor detail becomes important as will be discussed later in Section 8.

These considerations will be generalized for a particle obeying fractional dynamics. It can then be extended to the case of a particle in a potential field by including the potential term and making appropriate changes as will be described later.

3. Fractional Dynamics of a Free Particle

The first step is to generalize the equations of motion, (1)–(3), by replacing the integrals and derivatives by appropriate fractional integrals and derivatives. For physical problems with well-definable initial conditions the accepted practice is to employ the Caputo fractional derivatives [5]. The Caputo derivative of order 𝛽 is defined by [3]: ^𝐶₀𝐷^𝛽_𝑡𝑓(𝑡)=1Γ𝑛−𝛽^𝑡₀𝑓^𝑛(𝜏)𝑑𝜏(𝑡−𝜏)^{𝛽+1−𝑛}𝑛−1<𝛽<𝑛. (25) In the limit 𝛽→𝑛, the Caputo derivative becomes the ordinary 𝑛th derivative of the function.

The fractional integral of order 𝛽 is defined by ₀𝐼^𝛽_𝑡𝑓(𝑡)=1Γ𝛽^𝑡₀𝑓(𝜏)(𝑡−𝜏)^𝛽−1𝑑𝜏. (26) In generalizing the equations of motion, the second-order time derivative in Newton’s law is replaced by a Caputo derivative of order 𝛼, and the first-order derivative is replaced by a Caputo derivative of order (𝛼/2) [47]. Then (2) and (3) become 𝑥=𝑥₀+1𝑚_𝑓Γ(𝛼/2)^𝑡₀𝑝_𝑓(𝜏)(𝑡−𝜏)^𝛼/2−1𝑝𝑑𝜏,_𝑓=𝑝_𝑓0+1Γ(𝛼/2)^𝑡₀𝐹(𝜏)(𝑡−𝜏)^𝛼/2−1𝑑𝜏, (27) where as usual 𝑥₀, 𝑝_𝑓0 refer to the initial position and initial value of the “fractional momentum” 𝑝_𝑓, respectively. It is to be noted that the variables 𝑥, 𝑡 are still the space and time variables and have the dimensions of length and time [𝐿] and [𝑇], respectively. However, the dynamical quantities obtained by the operation of fractional derivation have different dimensions; for example, “fractional velocity” with the notation ^𝐶₀𝐷_𝑡^𝛼/2𝑥=𝑥̇^𝛼/2 would have the units [𝐿/𝑇^𝛼/2]. The dimension for the parameter “𝑚_𝑓” in the fractional momentum, 𝑝_𝑓=𝑚_𝑓𝑥̇^𝛼/2, is no longer just [𝑀] but has to be chosen [47] so that the fractional quantity 𝑝_𝑓²/2𝑚_𝑓 has the dimensions of energy [𝑀𝐿²/𝑇²]. Thus the dimension of the parameter “𝑚_𝑓” is [𝑀𝑇^2−𝛼] and the fractional momentum has the dimensions [𝑀𝐿/𝑇^2−𝛼/2]. The Lagrangian in (4) when generalized has the dimensions of energy. Of course, all quantities regain the standard dimensions in the limit 𝛼→2.

4. Time Fractional Schrodinger Equation for a Free Particle

There are two possible generalizations of the action integral in (5) used in fractional dynamics [11]: 𝑆_I=^𝑡₀𝐿𝑥,𝑥̇^𝛼/2𝑆,𝑡𝑑𝑡,_II=1Γ(𝛼/2)^𝑡₀𝐿𝑥,𝑥̇^𝛼/2,𝜏(𝑡−𝜏)^𝛼/2−1𝑑𝜏. (28)(29) Since the Lagrangian has the dimensions of energy, the dimensions of action defined in (28) and (29), 𝑆_I and 𝑆_II, are different.

The dimensions of 𝑆_I are the same as that of the regular action, namely, [𝑀𝐿²/𝑇], but that of 𝑆_II is [𝑀𝐿²/𝑇^2−𝛼/2]. In the Newtonian limit 𝛼→2, 𝑆_II→ regular action. These dimensional considerations have to be kept in mind in generalizing the Feynman method. In particular, if the choice from (29) is made, then a “fractional Planck constant” ℏ_𝑓 with appropriate dimensions must be introduced in order to render the argument of the exponential in (7) dimensionless.

For a “free particle” (27) yields 𝑥=𝑥₀+𝑝_𝑓0𝑡^𝛼/2𝑚_𝑓Γ(),𝑝1+𝛼/2_𝑓=𝑝_𝑓0. (30)(31) After carrying out the time slicing operation as in (11) and making the same approximation the evolution of the wave function in an infinitesimal interval of time 𝜀 can now be obtained. Equation (30) yields 𝑝_𝑓0=𝑥−𝑥₀𝑚_𝑓Γ(1+𝛼/2)𝜀^𝛼/2 (32) and hence 𝑚𝐿=^𝑓Γ(1+𝛼/2)𝑥−𝑥₀²2𝜀^𝛼. (33) But for action, there are two choices: 𝑆_I=𝑚_𝑓Γ²(1+𝛼/2)𝑥−𝑥₀²2𝜀^𝛼−1,𝑆_II=𝑚_𝑓Γ(1+𝛼/2)𝑥−𝑥₀²2𝜀^𝛼/2. (34)(35)

Making the appropriate changes, the equations for the evolution of the wave function in the two cases are 𝜓_I()=𝑥,𝑡+𝜀_∞−∞1𝐴_Iexp𝑖𝑚_𝑓Γ²(1+𝛼/2)𝑥−𝑥₀²2ℏ𝜀^𝛼−1×𝜓_I𝑥₀,𝑡𝑑𝑥₀,𝜓_II()=𝑥,𝑡+𝜀_∞−∞1𝐴_IIexp𝑖𝑚_𝑓Γ(1+𝛼/2)(𝑥−𝑥₀)²2ℏ_𝑓𝜀^𝛼/2×𝜓_II𝑥₀,𝑡𝑑𝑥₀. (36)(37) Changing the variable in the integrals from 𝑥₀ to 𝜂=𝑥−𝑥₀ as before and introducing two constants 𝑎_I=𝑚_𝑓Γ²(1+𝛼/2)2𝑖ℏ,𝑎_II=𝑚_𝑓Γ(1+𝛼/2)2𝑖ℏ_𝑓. (38) Equations (36) and (37) can be written as 𝜓_I()=𝑥,𝑡+𝜀_∞−∞1𝐴_I−𝑎exp^I𝜂²𝜀^𝛼−1𝜓_I𝜓𝑥+𝜂,𝑡𝑑𝜂,_II()=𝑥,𝑡+𝜀_∞−∞1𝐴_II−𝑎exp^II𝜂²𝜀^𝛼/2𝜓_II𝑥+𝜂,𝑡𝑑𝜂. (39)(40) The left-hand sides of (39) and (40) can be expanded in fractional Taylor’s series [48] in time, with fractional derivative of order 𝛾, and keeping only the lowest-order term in 𝛾 yields 𝜓_I,II(𝑥,𝑡+𝜀)=𝜓_I,II()+𝑥,𝑡_𝐶0𝐷^𝛾_𝑡𝜓_I,II()𝜀𝑥,𝑡^𝛾Γ𝛾+1+⋯. (41) In the right-hand side a Taylor expansion with terms up to second order in 𝜂 with respect to space can be carried out and the two cases will be considered separately. Thus, the equation for 𝜓_I(𝑥,𝑡) becomes 𝜓_I(𝑥,𝑡)+^𝐶₀𝐷^𝛾_𝑡𝜓_I()𝜀𝑥,𝑡^𝛾Γ=𝛾+1_∞−∞1𝐴_I−𝑎exp^I𝜂²𝜀^𝛼−1×𝜓_I(𝑥,𝑡)+𝜂𝜕𝜓_I+𝜂𝜕𝑥2²𝜕²𝜓_I𝜕𝑥²𝑑𝜂. (42) On evaluating the integrals on the right-hand side of (42), the middle term with the first power of 𝜂 vanishes. Equating the leading terms on both sides of (42) yields

𝜓_I(𝑥,𝑡)=(1/𝐴_I)𝜋𝜀^𝛼−1/𝑎_I𝜓_I(𝑥,𝑡), requiring that 𝐴_I=𝜋𝜀^𝛼−1/𝑎_I. Equation (42) reduces to 𝜓_I(𝑥,𝑡)+^𝐶₀𝐷^𝛾_𝑡𝜓_I()𝜀𝑥,𝑡^𝛾Γ=𝜓𝛾+1_I()+1𝑥,𝑡4𝜀^𝛼−1𝑎_I𝜕²𝜓_I(𝑥,𝑡)𝜕𝑥². (43) Equating the remaining terms requires that the powers of 𝜀 must be the same on both sides; that is, 𝛾=𝛼−1. Inserting the value of 𝑎_I and simplifying (43) yield ^𝐶₀𝐷_𝑡^𝛼−1𝜓_I(𝑥,𝑡)=𝑖ℏ2𝑚_𝑓Γ(𝛼)Γ²()𝜕1+𝛼/2²𝜓_I(𝑥,𝑡)𝜕𝑥². (44) Similarly, by expanding (40) the equation for 𝜓_II(𝑥,𝑡) becomes 𝜓_II(𝑥,𝑡)+^𝐶₀𝐷^𝛾_𝑡𝜓_II()𝜀𝑥,𝑡^𝛾Γ=𝛾+1_∞−∞1𝐴_II−𝑎exp^II𝜂²𝜀^𝛼/2×𝜓_II(𝑥,𝑡)+𝜂𝜕𝜓_II+𝜂𝜕𝑥2²𝜕²𝜓_II𝜕𝑥²𝑑𝜂. (45) On evaluating the integrals on the right-hand side of (45) the middle term vanishes. Equating the leading terms on both sides of (45) as before yields for the normalizing factor 𝐴_II=𝜋𝜀^𝛼/2/𝑎_II. Equation (45) reduces to 𝜓_II(𝑥,𝑡)+^𝐶₀𝐷^𝛾_𝑡𝜓_II()𝜀𝑥,𝑡^𝛾Γ=𝜓𝛾+1_II()+1𝑥,𝑡4𝜀^𝛼/2𝑎_II𝜕²𝜓_II(𝑥,𝑡)𝜕𝑥². (46) In the remaining terms, the powers of 𝜀 must be the same. This requires 𝛾=𝛼/2. Inserting the value of 𝑎_II and simplifying (46) yield the equation for the wave function 𝜓_II(𝑥,𝑡): ^𝐶₀𝐷_𝑡^𝛼/2𝜓_II(𝑥,𝑡)=𝑖ℏ_𝑓2𝑚_𝑓𝜕²𝜓_II(𝑥,𝑡)𝜕𝑥². (47) This completes the derivation based on the Feynman path integral method and (44) and (47) constitute the time fractional Schrodinger equations for a free particle corresponding to two ways of defining the action integral for a fractional dynamical system. In the limit 𝛼→2, the fractional dynamical system goes over to the regular Newtonian system and in this case both (44) and (47) reduce to 𝜕𝜓(𝑥,𝑡)=𝜕𝑡𝑖ℏ𝜕2𝑚²𝜓(𝑥,𝑡)𝜕𝑥² (48) given earlier, thus recovering the standard Schrodinger equation for a free quantum particle [49].

It should be noted that since 𝛼≤2, with 𝛼=2 being the limiting case, the order of the time fractional derivative is ≤1 in both (44) and (47) and cannot exceed 1. This means the time fractional Schrodinger equation as derived from the path integral method always corresponds to the “subdiffusion” case in contrast to the case where the TFSE is obtained by a simple replacement of the first-order time derivative by a fractional order derivative [19]. Furthermore, the order 𝛼/2 of the fractional derivative corresponds to the first-order regular derivative as has been used above in Section 3. Thus it appears that the second method of defining the action leading to (47) is the natural way to generalize to TFSE. In appearance also it is as if the equation has been obtained by replacing all quantities in the Schrodinger equation by an equivalent fractional quantity, except the space derivative. Furthermore, (44) becomes a fractional order integro-differential equation when 𝛼<1 and no longer just a fractional order differential equation. Because of these reasons, it is considered the method of choice to use (47) as the TFSE derived from the path integral method and no further reference will be made to (44).

The coefficient of the space derivative term in (47) has the dimension [𝐿²/𝑇^𝛼/2], corresponding to the fractional diffusion coefficient. Thus (47) can be considered a time fractional diffusion equation with an imaginary fractional diffusion coefficient, just as (48), the regular Schrodinger equation, can be considered to be a diffusion equation with an imaginary diffusion coefficient. Thus all the mathematical machinery of time fractional diffusion theory [5, 50–59] can be imported advantageously.

Although it is possible to introduce additional parameters and cast (47) in a dimensionless form, it has not been done here. However, for convenience, the subscript II is dropped from the wave function; a simplified notation for the Caputo derivative, 𝜕^𝛽_𝑡𝜓, with 𝛽=𝛼/2 will be used. Naturally, 0<𝛽≤1, and 𝛽→1 yields the regular first-order time derivative, denoted by 𝜕_𝑡. After incorporating these changes and defining a new constant 𝐷_𝑓=ℏ_𝑓/2𝑚_𝑓, (47) becomes 𝜕^𝛽_𝑡𝜓(𝑥,𝑡)=𝑖𝐷_𝑓𝜕²𝜓(𝑥,𝑡)𝜕𝑥²0<𝛽≤1. (49)

Equation (49) can be solved by a combination of the Fourier and Laplace transform methods [52–54].

5. Probability Current and the Continuity Equation

The probability density 𝜌 is defined by 𝜌=𝜓^∗𝜓=𝜓𝜓^∗. The complex conjugate wave function satisfies 𝜕^𝛽_𝑡𝜓^∗(𝑥,𝑡)=−𝑖𝐷_𝑓𝜕²𝜓^∗(𝑥,𝑡)𝜕𝑥². (50) There is an identity satisfied by the Caputo derivative [5] 𝜕_𝑡^1−𝛽𝜕^𝛽_𝑡𝑓(𝑡)=𝜕_𝑡𝑓(𝑡), (51) where the right-hand side represents the regular first-order derivative. This identity can be used in studying the time derivative of the probability density, given by 𝜕_𝑡𝜌=𝜕_𝑡𝜓𝜓^∗=𝜕_𝑡𝜓𝜓^∗+𝜓𝜕_𝑡𝜓^∗. (52) Inserting from (51) gives 𝜕_𝑡𝜌=𝜕_𝑡^1−𝛽𝜕^𝛽_𝑡𝜓𝜓^∗+𝜓𝜕_𝑡^1−𝛽𝜕^𝛽_𝑡𝜓^∗. (53) Substituting from (49) and (50), (53) yields 𝜕_𝑡𝜌=𝜕_𝑡^1−𝛽𝑖𝐷_𝑓𝜕²𝜓(𝑥,𝑡)𝜕𝑥²𝜓^∗+𝜓𝜕_𝑡^1−𝛽−𝑖𝐷_𝑓𝜕²𝜓^∗(𝑥,𝑡)𝜕𝑥². (54) Equation (54) can be rewritten after factoring out the constant and interchanging the order of space and time derivatives as 𝜕_𝑡𝜌=𝑖𝐷_𝑓𝜕²𝜕𝑥²𝜕_𝑡^1−𝛽𝜓𝜓^∗𝜕−𝜓²𝜕𝑥²𝜕_𝑡^1−𝛽𝜓^∗. (55) Introducing a new function 𝜓=𝜕_𝑡^1−𝛽𝜓, (55) can be written as 𝜕_𝑡𝜌=𝑖𝐷_𝑓𝜕²𝜕𝑥²𝜓𝜓^∗𝜕−𝜓²𝜕𝑥²𝜓^∗. (56) Defining a probability current density given by 𝐽_𝑥=−𝑖𝐷_𝑓𝜓∗𝜕𝜓𝜕𝜓𝜕𝑥−𝜓^∗𝜕𝑥, (57) equation (56) can be written as 𝜕_𝑡𝜌+𝜕𝐽_𝑥𝜕𝑥=0. (58) This is the time fractional version of the continuity equation. In the limit 𝛽→1, 𝜓→𝜓  𝜓^∗→𝜓^∗ and (58) reproduce the continuity equation of standard quantum mechanics [49]. It may be noted that (58) differs from Naber’s result ((24) in [19]) and will be discussed later.

6. Solution for TFSE for a Free Particle

The solution for the TFSE for a free particle under the conditions

𝜓(𝑥,0)=𝜓₀(𝑥); 𝜓(𝑥,𝑡)→0, |𝑥|→∞, 𝑡>0 is available in the literature but considered here for purposes of obtaining the Green function for the TFSE.

By applying the combined Fourier and Laplace transforms defined by 𝜓()=1𝑘,𝑠√2𝜋^∞−∞𝑒^{−𝑖𝑘𝑥∞}₀𝑒^−𝑠𝑡𝜓(𝑥,𝑡)𝑑𝑡𝑑𝑥 (59) equation (49) reduces to 𝑠^𝛽𝜓(𝑘,𝑠)−𝑠^𝛽−1𝜓(𝑘,0)=−𝑖𝐷_𝑓𝑘²𝜓(𝑘,𝑠) (60) resulting in 𝜓()=𝜓()𝑠𝑘,𝑠𝑘,0^𝛽−1𝑠^𝛽+𝑖𝐷_𝑓𝑘². (61)

Applying the inverse Laplace transform, (61) yields 𝜓(𝑘,𝑡)=𝜓(𝑘,0)𝐸_𝛽,1−𝑖𝐷_𝑓𝑘²𝑡^𝛽 (62) in terms of the Mittag-Leffler function defined by a series or by the inverse Lapalce transform [3]: 𝐸_𝛼,𝛽(𝑧)=^∞_𝑛=0𝑧^𝑛Γ𝛼𝑛+𝛽=𝐿⁻¹𝑠^𝛼−𝛽𝑠^𝛼.−𝑧 (63)

The Green function solution can be written as 𝜓(𝑥,𝑡)=^∞_−∞𝜓(𝑥−𝜉)𝐺_𝛽(𝜉,𝑡)𝑑𝜉, (64) where the Green function is given by the inverse Fourier transform of the Mittag-Leffler function 𝐺_𝛽()=1𝑥,𝑡√2𝜋^∞−∞𝑒^𝑖𝑘𝑥𝐸_𝛽,1−𝑖𝐷_𝑓𝑘²𝑡^𝛽𝑑𝑘, (65) where the Mittag-Leffler function has been defined [3] before in (63).

The Fourier inversion in (65) can be carried out [52–54] using the property that the Mittag-Leffler function is related through the Laplace integral to another special function of the Wright type denoted by 𝑀_𝛽(𝑧)==𝑊−𝑧;−𝛽,1−𝛽_∞𝑛=0(−𝑧)^𝑛𝑛!Γ−𝛽𝑛+1−𝛽0<𝛽<1, (66) where the Wright function is defined by [54] 𝑊𝑧;𝜆,𝜇=^∞_𝑛=0𝑧^𝑛𝑛!Γ𝜆𝑛+𝜇𝜆>−1,𝜇∈𝐶. (67) The Green function in (65) is then given by 𝐺_𝛽()=1𝑥,𝑡21𝑖𝐷_𝑓𝑡^𝛽/2𝑀_𝛽/2|𝑥|𝑖𝐷_𝑓𝑡^𝛽/2. (68) In the context of fractional diffusion, the function 𝑀_𝛽/2(𝑧) belongs to the Wright type of probability densities characterized by the similarity variable 𝑧=|𝑥|/𝐷_𝑓𝑡^𝛽/2, where 𝐷_𝑓 is the fractional diffusion coefficient, with ∫^∞₀𝑀_𝛽/2(𝑧)𝑑𝑧=1. Furthermore, the probability densities are non-Markovian and exhibit a variance consistent with slow anomalous diffusion [49–51], 𝜎²_𝛽(𝑡)=(2/Γ(𝛽+1))𝐷_𝑓𝑡^𝛽.

In the limit of 𝛽→1, the probability density function goes over to the Gaussian 𝑀_1/2(𝑧)=1√𝜋𝑒^−𝑧2/4 (69) corresponding to regular diffusion. A plot of the reduced probability density function is given in Figure 1.

Figure 1 

Reduced probability density function.

The Green function for regular diffusion describes a probability density, whereas the corresponding Green function for the Schrodinger equation is the propagator, which describes the probability amplitude for the particle to propagate from 𝑥_𝑎 at 𝑡_𝑎 to 𝑥_𝑏 at 𝑡_𝑏. In exactly the same way, the Green function for time fractional diffusion describes a probability density, whereas the Green function in (68) is the fractional propagator and gives the probability amplitude. Of particular interest is the Fourier component of the wave function in (62) in connection with the total probability as 𝑡→∞; the case discussed by Naber [19] and will be discussed in the next section.

7. Comments about Some Results in Naber’s Work [19]

This section draws attention to some errors in Naber’s paper which have gone unnoticed and have been reproduced repeatedly. Naber’s equations will be referred to by their number and a prefix N. Naber explicitly states that he uses the Caputo derivative, which has been defined earlier in (25) in this paper. Although Naber does not give the explicit definition of the Caputo fractional derivative in his paper [19], it can be inferred from (NA.3) given in the appendix to his paper.

(a) However, Naber writes in (N16) for a Caputo derivative of order (1−𝜈), reproduced here for convenience in (70): 𝐷_𝑡^1−𝜈𝜓()=1𝑡,𝑥Γ(1−𝜈)×^𝑡₀𝑑𝜓(𝜏,𝑥)𝑑𝜏𝑑𝜏(𝑡−𝜏)^𝜈(0<𝜈<1). (70) This is incorrect, as can be checked easily by taking the limit 𝜈→1. The left-hand side →𝜓(𝑡,𝑥) as the derivative of zero order, but the right-hand side →0 because of the Γ function in the denominator. The correct form of equation is 𝐷_𝑡^1−𝜈𝜓()=1𝑡,𝑥Γ(𝜈)^𝑡₀𝑑𝜓(𝜏,𝑥)𝑑𝜏𝑑𝜏(𝑡−𝜏)^1−𝜈. (71) As a consequence, the weight factor in (N18) should be (𝑡−𝜏)^1−𝜈 and not (𝑡−𝜏)^−𝜈, which Naber uses to give a physical significance to the entity he has introduced. The correct weight factor in (N18) would →1 in the limit 𝜈→1.

(b) In (N11), Naber gives an identity satisfied supposedly by fractional Caputo derivatives of order less than 1, reproduced here for convenience in (72). 𝐷_𝑡^1−𝜈𝐷^𝜈_𝑡𝑦(𝑡)=𝑑𝑦−𝐷𝑑𝑡^𝜈𝑡𝑦(𝑡)_𝑡=0𝑡^1−𝜈Γ(𝜈). (72) This identity is not correct and cannot be found anywhere. The identity satisfied by the Caputo derivatives is given in [5] and reproduced later in the current notation for 𝜈<1𝐷_𝑡^−𝜈𝐷^𝜈_𝑡𝑦(𝑡)=𝑦(𝑡)−𝑦(0). (73) This yields 𝐷_𝑡^1−𝜈𝐷^𝜈_𝑡𝑦(𝑡)=𝑑𝑦/𝑑𝑡=𝜕_𝑡𝑦(𝑡) and has been used earlier in (52) in this paper.

The incorrect identity has been used by Naber to derive an equation for the probability current, which is obviously incorrect. Unfortunately, the incorrect identity, as given by Naber, in (72) has been repeatedly used in the literature. The correct equation for the probability current has been given in this paper in (59), which reduces to the standard continuity equation for the probability current in regular quantum mechanics [49].

(c) This point concerns the separation of the Mittag-Leffler function with an imaginary argument into an oscillatory part and a part which decays exponentially with time. There is nothing wrong with the derivation itself as given by Naber, and the function under discussion is the Fourier component of the free particle wave function, in his notation Ψ=Ψ₀𝐸_𝜈𝜔(−𝑖𝑡)^𝜈, (74) where 𝐸_𝜈(𝑧) in Naber’s notation corresponds to the Mittag-Leffler function 𝐸_𝜈,1(𝑧) defined in (63).

The Mittag-Leffler function with the complex argument has been separated into an oscillatory part and a part based on the evaluation of the inverse Laplace transform 𝐴(𝑡)=12𝜋𝑖^𝛾+𝑖∞_{𝛾−𝑖∞}𝑒^𝑠𝑡𝑠^𝜈−1𝐴₀𝑑𝑠𝑠^𝜈−𝜎𝑖^𝜈 (75) along a Hankel contour and considering the contribution from the residue of the pole 𝑠₀=𝜎^1/𝜈𝑖 together with the contribution from the integral along the two strips on either side of the branch cut, which is a standard procedure. Naber finally gives the solution as Ψ=Ψ₀𝑒^{−𝑖𝜔1/𝜈𝑡}𝜈−𝐹_𝜈𝜔(−𝑖)^𝜈,𝑡. (76) He argues that in the limit of 𝑡→∞ the total probability arises basically from the first term and is equal to 1/𝜈², assuming that the wave function was initially normalized. Since 𝜈<1, the total probability is >1, a result difficult to understand physically. However, it is shown later that the solution derived in this paper yields a probability that →0 in the limit 𝑡→∞.

The solution corresponding to (74) is given by (62) in this paper. The separation into two parts corresponding to (76) can be carried out by considering the inverse Laplace transform of (61): 𝜓()=1𝑘,𝑡2𝜋𝑖^𝛾+𝑖∞_{𝛾−𝑖∞}𝑒^𝑠𝑡𝑠^𝛽−1𝜓(𝑘,0)𝑑𝑠𝑠^𝛽+𝑖𝐷_𝑓𝑘². (77) As usual, the Bromwich contour is replaced by the Hankel contour. The two contributions arise from (i) the residue at the pole at 𝑠₀=(𝐷_𝑓𝑘²)^1/𝛽(−𝑖)^1/𝛽 and (ii) the integral along the two strips from 0 to −∞ introduced along the branch cut. The latter yields a contribution which decays in time just as the second term in Naber’s equation (76). The contribution from the residue is given by Residue𝜓()=𝑒𝑘,0𝛽^{(−𝑖)1/𝛽(𝐷𝑓𝑘2)1/𝛽𝑡} (78) and corresponds to the first term on the right-hand side in (76). However, (−𝑖)^1/𝛽=𝑒^{−𝑖𝜋/21/𝛽}=𝑒^{−𝑖𝜋/2𝛽}=𝜋cos𝜋2𝛽−𝑖sin.2𝛽 (79) Therefore, (78) yields 𝜓(𝑘,𝑡)𝜓()=𝑘,0Residue𝜓(𝑘,0)=𝑒^{𝑡(𝐷𝑓𝑘2)1/𝛽cos(𝜋/2𝛽)}𝑒^{−𝑖𝑡(𝐷𝑓𝑘2)1/𝛽sin(𝜋/2𝛽).} (80) The right-hand side of (80) has an amplitude term and an oscillatory factor. Because 𝛽<1, cos(𝜋/2𝛽) is negative, the amplitude factor decays in time exponentially and →0 in the limit 𝑡→∞. Thus the contribution from the residue due to the pole also →0. Thus the entire wave function →0 in the limit of 𝑡→∞, so does the probability density. Hence the total probability also →0, in contrast to Naber’s result.

The question naturally arises why there is this difference in the two results, both of which are concerned with the Mittag-Leffler functions with complex arguments. The reason appears to be that in Naber’s case, the pole occurs at 𝑠₀=𝜎^1/𝜈𝑖, whereas in the present paper, the pole occurs (using the same notation as Naber) at 𝑠₀=𝜎^1/𝜈(−𝑖)^1/𝜈. The simple 𝑖 in Naber’s case leads to the purely oscillatory solution and hence to the result of the probability being greater than unity. In our paper the imaginary number raised to the fractional power leads to the exponentially decaying solution and hence leads to the correct limit when 𝑡→∞, namely, zero total probability. Naber’s result is a direct consequence of his choice to raise the power of the imaginary number 𝑖 to the fractional power, so as to incorporate a Wick rotation. However, as had been indicated earlier, this imaginary number cannot be arbitrarily altered as it is connected with the phase factor 𝑖𝑆/ℏ in the Feynman propagator. If it is necessary to include a Wick rotation, the power of 𝑖 should be changed to 𝜈+1 in Naber’s equation (N9) instead of just 𝜈. If this is done, the total probability would properly →0 in the limit 𝑡→∞. If this is done, then the TFSE can be interpreted as the analytic continuation of the fractional diffusion equation, just as the regular Schrodinger equation can be considered as the analytic continuation of the regular diffusion equation.

8. TFSE for a Particle in a Potential Field

So far attention has been focused on a free particle. These considerations can be easily extended to the case of a particle moving in a potential field by incorporating a potential term 𝑉(𝑥,𝑡) in the Lagrangian 𝐿=𝑇−𝑉. This will necessitate incorporating an additional term −𝑉(𝑥,𝑡)𝜀^𝛼/2/Γ(1+𝛼/2) in the right-hand side of (35) and an additional factor −(𝑖/ℎ_𝑓)(𝑉(𝑥,𝑡)𝜀^𝛼/2/Γ(1+𝛼/2)) in the exponential in (37). This results in changing (45) into 𝜓_II(𝑥,𝑡)+^𝐶₀𝐷^𝛽_𝑡𝜓_II()𝜀𝑥,𝑡^𝛽Γ=𝛽+1_∞−∞1𝐴_II−𝑎exp^II𝜂²𝜀^𝛼/2𝑖1−ℎ_𝑓𝑉(𝑥,𝑡)𝜀^𝛼/2Γ()×𝜓1+𝛼/2_II(𝑥,𝑡)+𝜂𝜕𝜓_II+𝜂𝜕𝑥2²𝜕²𝜓_II𝜕𝑥²𝑑𝜂. (81)

Equation (47) then becomes ^𝐶₀𝐷_𝑡^𝛼/2𝜓_II(𝑥,𝑡)=𝑖ℏ_𝑓2𝑚_𝑓𝜕²𝜓_II(𝑥,𝑡)𝜕𝑥²−𝑖ℎ_𝑓𝑉(𝑥,𝑡)𝜓(𝑥,𝑡). (82) Multiplying both sides by −ℏ_𝑓/𝑖 results in the FTSE for a particle in a potential field as 𝑖ℏ_𝑓𝜕^𝛽_𝑡𝜓()ℏ𝑥,𝑡=−^2𝑓2𝑚_𝑓𝜕²𝜓(𝑥,𝑡)𝜕𝑥²−𝑉(𝑥,𝑡)𝜓(𝑥,𝑡) (83) which reduces to the standard Schrodinger equation of quantum mechanics: 𝑖ℏ𝜕𝜓ℏ𝜕𝑡=−²𝜕2𝑚²𝜓𝜕𝑥²+𝑉𝜓. (84) One final remark needs to be made concerning the use of the Planck units for casting the TFSE in terms of dimensionless quantities. After the derivation, the FTSE in (83) can be cast in dimensionless quantities; however, the appropriate fractional Planck units must be defined to take care of the fractional quantities involved in (83).

9. Discussion and Conclusions

The TFSE has been derived using the Feynman path integral technique for a nonrelativistic particle. As expected the TFSE looks like a time fractional diffusion equation with an imaginary fractional diffusion constant but pertains to the realm of subdiffusion only, in contrast to Naber’s generalization which includes superdiffusion as well. This is understandable because the case considered in this paper pertains to the nonrelativistic case. Relativistic considerations would have to be included for the superdiffusion case, which would lead to the Klein-Gordon equation in the integer order limit. In adapting the Feynman method, it is shown that it is preferable to introduce the action 𝑆 as a fractional time integral of the Lagrangian and that it is necessary to introduce a “fractional Planck constant.” In the limit of integer order, the regular action 𝑆, the regular Schrodinger equation, and the regular Planck constant are all recovered. It may be of interest to note that there is a fractional Planck constant implied in Naber’s work also, although it is not explicitly so stated. His equations are rendered nondimensional by using the Planck units of mass, length, and time and then generalized to fractional derivatives after including a Wick rotation. This implies a change of variable of time 𝑡→𝑖𝑡 so that the imaginary number is raised to the same power as the order of the fractional time derivative involved. However, the Planck units may not be the appropriate quantities as the equations involve quantities of fractional dimension and the equations must be made nondimensional after the generalization to fractional derivatives and not before. A fractional Planck constant does show up in Naber’s treatment as well, as the ratio of masses 𝑚/𝑀_𝑝 [31, 32].

A number of errors in Naber’s work have been corrected. The correct continuity equation for the probability current is derived and the Green function solution for a free particle is given. The Green function is given in terms of a special type of function, the M-Wright function, which is used extensively in studies of time fractional subdiffusion studies. In the context of time fractional diffusion, the M-Wright function is a probability density function in time, which is non-Markovian and goes over to the Gaussian in the nonfractional limiting case. In the context of TFSE, the M-Wright function gives the propagator, which is a probability amplitude. Probability considerations are accounted for by the usual process of squaring of the amplitude. In particular, it is shown that Naber’s result that the total probability is greater than unity in the long-time limit is a spurious result arising out of the operation of arbitrarily raising of the imaginary number 𝑖 to the power of the same degree as the fractional time derivative invoking a Wick rotation. It is shown that such arbitrary change of the imaginary number cannot be carried out as the imaginary number 𝑖 is connected with the phase of the action 𝑆 in the path integral contribution. However, if we desired to consider a Wick rotation, it should be included as an additional increase of the index of the power of the imaginary number. The TFSE for a particle moving in a potential field is also derived. Furthermore, it is suggested that even in studies of fractional classical mechanics, such as those using variational methods, the action integral be expressed as a time fractional integral of the Lagrangian. Further extensions including the solutions to particle subject to different potentials are underway.