Stability of Perfectly Matched Layers for Time Fractional Schrödinger Equation

Tingting Zhang^*, Xiangkun Li

School of Mathematics and Statistics, Shandong Normal University, Jinan, China.

DOI: 10.4236/eng.2023.151001                

Abstract

It is an important issue to numerically solve the time fractional Schrödinger equation on unbounded domains, which models the dynamics of optical solitons propagating via optical fibers. The perfectly matched layer approach is applied to truncate the unbounded physical domain, and obtain an initial boundary value problem on a bounded computational domain, which can be efficiently solved by the finite difference method. The stability of the reduced initial boundary value problem is rigorously analyzed. Some numerical results are presented to illustrate the accuracy and feasibility of the perfectly matched layer approach. According to these examples, the absorption parameters and the width of the absorption layer will affect the absorption effect. The larger the absorption width, the better the absorption effect. There is an optimal absorption parameter, the absorption effect is the best.

Keywords

Time Fractional Schrödinger Equation, Perfectly Matched Layer, Stability

1. Introduction

As one of the most important model in modern science, the Schrödinger equation has a variety of applications in many fields, such as optics, Bose-Einstein condensation, quantum mechanics and molecular physics. The Schrödinger equation can be derived by considering Gaussian probability distribution [1]. It makes sense to consider the non-Gaussian distribution, and Laskin [2] [3] obtained the space fractional Schrödinger equation. The diffusion of particles is still Markovian for the space fractional Schrödinger equation, and it is no need to consider the memory of the equation. However, the memory and hereditary properties must be considered under the non-Markovian evolution. Naber firstly proposed the time fractional Schrödinger equation based on the generation of the time fractional diffusion equation [4], and studied the properties of the time fractional Schrödinger equation. In addition, Sjögreen proposed the perfectly matched layer method of Maxwell’s equations [5]. Singer and Turkel studied the Helmholtz equation based on the perfectly matched layer method [6]. Karim and the coauthors combined perfectly matched layer with finite element method [7] to study the numerical modeling of acoustic devices in piezoelectric. As a fundamental model of fractional quantum mechanics [8], the fractional Schrödinger equation [9] [10] has attracted massive attention from both physicists and mathematicians; see the references [11] [12] [13] and references therein. Although the approximate and analytical solutions have been studied for some initial conditions or nonlinear terms, the expression of the solution is very complex [14] [15]. Thus, how to design the efficient method to solve the fractional Schrödinger equation on unbounded domain is an important issue.

This paper aims to numerically solve the following time fractional Schrödinger equation (TFSE) on unbounded domain

$i{}_0{}^{C}D{}_t^{\alpha }\psi \left(x,t\right)=-{\psi }_{xx}\left(x,t\right),\mathrm{<mtext>\; \hspace{0.17em}\; </mtext>\; <mtext>\; \hspace{0.17em}\; </mtext>}\left(x,t\right)\in ℝ\times \left(0,T\right],$ (1)

$\psi \left(x,0\right)={\psi }_0\left(x\right),\mathrm{<mtext>\; \hspace{0.17em}\; </mtext>}x\in ℝ,$ (2)

$\psi \left(x,t\right)\to 0,\mathrm{<mtext>\; \hspace{0.17em}\; </mtext>}\left|x\right|\to \infty ,$ (3)

where $i=\sqrt{-1}$, $\psi \left(x,t\right)$ is the complex-valued wave function, ${}_0{}^{C}D{}_t^{\alpha }$ is the Caputo fractional derivative operator with $\alpha \in \left(\mathrm{0,1}\right)$ defined by

${}_0{}^{C}D{}_t^{\alpha }\psi \left(x,t\right)=\frac{1}{\Gamma \left(1-\alpha \right)}{\displaystyle {\int }_0^t}\frac{\partial \psi \left(x\mathrm{,}s\right)}{\partial s}\frac{1}{{\left(t-s\right)}^{\alpha }}\text{d}s$

with the usual Gamma function $\Gamma (\cdot )$. The function ${\psi }_0\left(x\right)$ is the initial condition with compactly supported in the domain of interest.

There are several approaches to efficiently solve the partial differential equations defined on unbounded domain, including infinite element method, boundary element method, artificial boundary method and absorbing layer approach. We restrict ourselves to the last strategy and numerically solve the TFSE by applying the perfectly matched layer (PML) technique, which was originally proposed by Berenger for electromagnetism [16] and has been successfully applied to solve a variety of partial differential equations [17] [18] [19], to study the behavior of quantum mechanical systems without having spurious reflections from waves traveling out of the interested domain. The key idea of PML approach is to surround the physical domain by an artificial unphysical layer to damp the waves entering the layer region without any reflections. Many works are presented to study the numerical solution of the Schrödinger-type equation by adopting the PML approach. Zheng applied the PML approach to solving the linear and nonlinear Schrödinger equation on unbounded domain in [20]. Nissen and Kreiss used a modal ansatz to derive the PML approach for the Schrödinger equation in [21], and demonstrated how to choose the optimal parameters of the PML. The PML approach for a system of two-dimensional coupled nonlinear Schrödinger equations with mixed derivatives is derived in [22], and the stability is proved for a specified absorption function. The recent theoretical and numerical developments concerning PML approach for solving the Schrödinger equation are presented in the review [23] and the references therein.

There are some works to study the numerical solution of the time fractional Schrödinger equation on unbounded domain. Li, Zhang and Antoine [24] [25] designed the artificial boundary conditions and fast algorithms to efficiently solve the TFSE in unbounded domains. Antoine and coauthors designed a series of PML for the time dependent space fractional partial differential equations [26] [27] [28]. Unfortunately, the study of the PML for TFSE is very rare. This paper aims to research this issue and the corresponding stability.

The rest of this paper is organized as follows. In Section 2, the general solution of the TFSE is obtained by applying the Laplace transformation, and the PML approach for TFSE is applied to obtain a reduced initial boundary value problem (IBVP) on a bounded domain. The stable property of IBVP with PML function is analyzed rigorously. The IBVP is discretized by applying the finite difference method in Section 3. In Section 4, some numerical results are presented to illustrate the accuracy and effectiveness of the our PML approach. Finally, the conclusions and research purposes are given in Section 5.

2. PML Approach for TFSE

2.1. PML Function

The perfectly matched layers with width d are introduced to divide the whole area into three parts, the interior domain ${\Omega }_i=\left[x_l\mathrm{,}{x}_r\right]$ with adjacent PML is reduced to the computational domain ${\Omega }_c=\left[x_l-d\mathrm{,}{x}_r+d\right]$, and exterior domain is ${\Omega }_e=ℝ\backslash {\Omega }_i$. In order to obtain the equation in the PML, which is modified through the idea of Nissen [21], consider the constraint on $\psi \left(x,t\right)$ on the exterior domain ${\Omega }_o=\left[x_r\mathrm{,}+\infty \right)$

$i\mathrm{<mtext>\; \hspace{0.17em}\; </mtext>}{}_0{}^{C}D{}_t^{\alpha }\psi \left(x\mathrm{,}t\right)=-{\psi }_{xx}\left(x\mathrm{,}t\right)\mathrm{,}x\in {\Omega }_o\times \left(\mathrm{0,}T\right]\mathrm{,}$ (4)

$\psi \left(x\mathrm{,0}\right)=\mathrm{0,}x\in {\Omega }_e\mathrm{,}$ (5)

$\psi \left(x\mathrm{,}t\right)\to \mathrm{0,}x\to \infty \mathrm{.}$ (6)

One notes that the Laplace transform of Caputo fractional derivative [29] satisfies

${}_0{}^{c}D{}_t^{\alpha }\left[\widehat{f\left(t\right)}\right]\left(s\right)=s^{\alpha }\hat{f}(s)-s^{\alpha -1}f\left(0\right).$

Applying the Laplace transform to Equation (4), we get

$-i{s}^{\alpha }\hat{\psi }\left(x\mathrm{,}s\right)={\hat{\psi }}_{xx}\left(x\mathrm{,}s\right)\mathrm{,}$

which has a general solution

$\hat{\psi }\left(x\mathrm{,}s\right)=c_1\left(s\right){\text{e}}^{-\sqrt[+]{-i{s}^{\alpha }}}+c_2\left(s\right){\text{e}}^{\sqrt[+]{-i{s}^{\alpha }}}\mathrm{.}$

Assume that $Re\left(\sqrt{-i{s}^{\alpha }}\right)\ge 0$ for $Re\left(s\right)\ge 0$. Since $\hat{\psi }\left(x\mathrm{,}s\right)\to 0$ when $x\to \infty$, we have

${\hat{\psi }}_{PML}\left(x\right)=c_1\exp \left\{-\sqrt[+]{-i{s}^{\alpha }}\left[\left(x-x_r\right)+{\text{e}}^{i\theta }{\displaystyle {\int }_{x_r}^x}\sigma \left(\omega \right)\text{d}\omega \right]\right\}\mathrm{,}$ (7)

where the absorption function $\sigma \left(\omega \right)$ is a real and non-negative function in $\omega$. To obtain decaying solution, the parameter $\theta$ is usually chosen as a constant

in the interval $\theta \in \left(\mathrm{0,}\frac{\pi }{2}\right)$. The PML function (7) satisfies

$i{s}^{\alpha }\hat{\psi }+\frac{1}{1+{\text{e}}^{i\theta }\sigma \left(x\right)}\frac{\partial }{\partial x}\left[\frac{1}{1+{\text{e}}^{i\theta }\sigma \left(x\right)}\frac{\partial \hat{\psi }}{\partial x}\right]=0.$ (8)

Applying the inverse Laplace transform to (8), we get the following initial boundary value problem in bounded domain

$i\mathrm{<mtext>\; \hspace{0.17em}\; </mtext>}{}_0{}^{C}D{}_t^{\alpha }\psi \left(x,t\right)+\gamma \left(x\right)\frac{\partial }{\partial x}\left[\gamma \left(x\right)\frac{\partial \psi }{\partial x}\right]=0,x\in {\Omega }_i\times \left(0,T\right],$ (9)

$\psi \left(x,0\right)={\psi }_0\left(x\right),\mathrm{<mtext>\; \hspace{0.17em}\; </mtext>\; <mtext>\; \hspace{0.17em}\; </mtext>}x\in {\Omega }_i,$ (10)

$\psi \left(x,t\right)=0,\mathrm{<mtext>\; \hspace{0.17em}\; </mtext>}\left(x,t\right)=\left\{x_l-d,x_r+d\right\}\times \left[0,T\right],$ (11)

where $\gamma \left(x\right)=\frac{1}{1+{\text{e}}^{i\theta }\sigma \left(x\right)}$. One can observe that $\sigma \left(x\right)=0$ for $x\in \left[x_l\mathrm{,}{x}_r\right]$ to

reduce the original Equation (1) in the bounded computational domain. Thus, Equation (9) can be solved both in the interior domain and in the layer by letting $\sigma \left(x\right)$ vanish in the interior.

2.2. Stability

Based on the idea of smooth exterior scaling (SES) [21], we present the stability of the reduced IBVP (9)-(11) in this section. In SES the ansatz $\psi \left(x\mathrm{,}t\right)=\kappa \left(x\right)\phi \left(x\mathrm{,}t\right)$ is introduced. Substituting the ansatz into (9), we have

$i\mathrm{<mtext>\; \hspace{0.17em}\; </mtext>}{}_0{}^{C}D{}_t^{\alpha }\phi \left(x\mathrm{,}t\right)=-\frac{1}{\upsilon \left(x\right)}\frac{{\partial }^2}{\partial x^2}\left[\frac{1}{\upsilon \left(x\right)}\phi \left(x\mathrm{,}t\right)\right]+V_{PML}\left(x\right)\phi \left(x\mathrm{,}t\right)\mathrm{,}$ (12)

where $\upsilon \left(x\right)=1+{\text{e}}^{i\theta }\sigma \left(x\right)$ and $V_{PML}\left(x\right)=\frac{3{\left[{\upsilon }^{\prime }\left(x\right)\right]}^2-2\upsilon \left(x\right){\upsilon }^{″}\left(x\right)}{4{\upsilon }^4\left(x\right)}$. The $\kappa \left(x\right)$ can be obtained as $\kappa \left(x\right)={\upsilon }^{-\frac{1}{2}}\left(x\right)$.

Lemma 1. [25] Let $u\left(t\right)$ be a complex valued function which is absolutely continuous on $\left[\mathrm{0,}T\right]$. Then, the following inequality holds

${}_0{}^{C}D{}_t^{\alpha }{\left|u\left(t\right)\right|}^2\le \bar{u}\left(t\right){}_0{}^{C}D{}_t^{\alpha }u\left(t\right)+u\left(t\right){}_0{}^{C}D{}_t^{\alpha }\bar{u}\left(t\right),\mathrm{<mtext>\; \hspace{0.17em}\; </mtext>\; 0}<\alpha <1,$ (13)

where $\bar{u}$ represents the complex conjugate function of u and ${\left|u\right|}^2=u\bar{u}$.

Lemma 2. [25] Let $0<\alpha <1$. We assume that $u\left(t\right)$ is a nonnegative absolutely continuous function that satisfies

${}_0{}^{C}D{}_t^{\alpha }u\left(t\right)\le c_{1}u\left(t\right)+c_2\left(t\right)\mathrm{,\; <mtext>\; \hspace{0.17em}\; </mtext>\; <mtext>\; \hspace{0.17em}\; </mtext>}t\in \left(\mathrm{0,}T\right]\mathrm{,}$ (14)

where $c_1>0$ and $c_2\left(t\right)$ is an integrable nonnegative function on $\left[\mathrm{0,}T\right]$. We have

$u\left(t\right)\le u\left(0\right)E_{\alpha \mathrm{,1}}\left(c_1{t}^{\alpha }\right)+\Gamma \left(\alpha \right)E_{\alpha \mathrm{,}\alpha }\left(c_1{t}^{\alpha }\right){}_{0}D{}_t^{-\alpha }{c}_2\left(t\right)\mathrm{,}$ (15)

where $E_{\alpha ,\mu }\left(z\right)={\displaystyle \underset{n=0}{\overset{\infty }{\sum }}}\frac{z^n}{\Gamma \left(n\alpha +\mu \right)}$ are the Mittag-Leffler functions and ${}_{0}D{}_t^{-\alpha }u\left(t\right)$ is the Riemann-Liouville integral given by ${}_{0}D{}_t^{-\alpha }u\left(t\right)=\frac{1}{\Gamma \left(\alpha \right)}{\displaystyle {\int }_0^t}\frac{u\left(\tau \right)}{{\left(t-\tau \right)}^{1-\alpha }}\text{d}\tau$.

Theorem 1. Let $\gamma \left(x\right)=\beta -i\eta$. Assume that $\beta \ge {\beta }_0>0$ , $\eta \ge 0$, $\left(\beta \mathrm{,}\eta \right)\in ℝ$ and that (12) has a smooth solution in any interval $0\le t\le T<\infty$. Then

${\Vert \phi \left(\cdot \mathrm{,}t\right)\Vert }^2\le K\left(T\right){\Vert \phi \left(\cdot \mathrm{,0}\right)\Vert }^2\mathrm{.}$ (16)

Proof. Since $\gamma \left(x\right)=\frac{1}{\upsilon \left(x\right)}$, the Equation (12) can be written as

$i\mathrm{<mtext>\; \hspace{0.17em}\; </mtext>}{}_0{}^{C}D{}_t^{\alpha }\phi \left(x\mathrm{,}t\right)=-\gamma \left(x\right){\left[\gamma \left(x\right)\phi \left(x\mathrm{,}t\right)\right]}_{xx}+V_{PML}\left(x\right)\phi \left(x\mathrm{,}t\right)\mathrm{,}$ (17)

where $V_{PML}\left(x\right)=\frac{2\gamma {\gamma }_{xx}-{\gamma }_x^2}{4}$.

By multiplying (17) by $\bar{\phi }\left(x\mathrm{,}t\right)$, taking the complex conjugate of (17), multiplying the result by $\phi \left(x\mathrm{,}t\right)$, and combining the two equations and integrating over ${\Omega }_i$, we get

$\begin{array}{l}{\displaystyle {\int }_{{\Omega }_i}}\left[\bar{\phi }\left(x\mathrm{,}t\right){}_0{}^{C}D{}_t^{\alpha }\phi \left(x\mathrm{,}t\right)+\phi \left(x\mathrm{,}t\right){}_0{}^{C}D{}_t^{\alpha }\bar{\phi }\left(x\mathrm{,}t\right)\right]\text{d}x\\ ={\displaystyle {\int }_{{\Omega }_i}}\left[i\left(\gamma \bar{\phi }\right){\left(\gamma \phi \right)}_{xx}-i\left(\bar{\gamma }\phi \right){\left(\bar{\gamma }\bar{\phi }\right)}_{xx}+i\left({\bar{V}}_{PML}-V_{PML}\right)\bar{\phi }\phi \right]\text{d}x\\ =I_1+I_2\mathrm{.}\end{array}$ (18)

$\begin{array}{c}{I}_1={\displaystyle {\int }_{{\Omega }_i}}\left[i\left(\gamma \bar{\phi }\right){\left(\gamma \phi \right)}_{xx}-i\left(\bar{\gamma }\phi \right){\left(\bar{\gamma }\bar{\phi }\right)}_{xx}\right]\text{d}x\\ ={\displaystyle {\int }_{{\Omega }_i}}\left[-i{\left(\bar{\phi }\left(\beta -i\eta \right)\right)}_x{\left(\phi \left(\beta -i\eta \right)\right)}_x+i{\left(\phi \left(\beta +i\eta \right)\right)}_x{\left(\bar{\phi }\left(\beta +i\eta \right)\right)}_x\right]\text{d}x\\ =-2{\displaystyle {\int }_{{\Omega }_i}}\left[{\left(\eta \bar{\phi }\right)}_x{\left(\beta \phi \right)}_x+{\left(\beta \bar{\phi }\right)}_x{\left(\eta \phi \right)}_x\right]\text{d}x.\end{array}$

Denote $\phi =z+iw$, we have

$\begin{array}{c}{I}_1=-2{\displaystyle {\int }_{{\Omega }_i}}\left[{\left(\eta \left(z-iw\right)\right)}_x{\left(\beta \left(z+iw\right)\right)}_x+{\left(\beta \left(z-iw\right)\right)}_x{\left(\eta \left(z+iw\right)\right)}_x\right]\text{d}x\\ =-4{\displaystyle {\int }_{{\Omega }_i}}\left[{\left(\beta z\right)}_x{\left(\eta z\right)}_x+{\left(\beta w\right)}_x{\left(\eta w\right)}_x\right]\text{d}x\\ =-4{\displaystyle {\int }_{{\Omega }_i}}\left[\beta \eta \left({\left|z_x\right|}^2+{\left|w_x\right|}^2\right)\right]\text{d}x-4{\displaystyle {\int }_{{\Omega }_i}}\left[\left({\eta }_x\beta +\eta {\beta }_x\right)\left(z{z}_x+w{w}_x\right)\right]\text{d}x\\ -4{\displaystyle {\int }_{{\Omega }_i}}\left[{\eta }_x{\beta }_x\left({\left|z\right|}^2+{\left|w\right|}^2\right)\right]\text{d}x\\ =-4{\displaystyle {\int }_{{\Omega }_i}}\left[\beta \eta \left({\left|z_x\right|}^2+{\left|w_x\right|}^2\right)\right]\text{d}x+4{\displaystyle {\int }_{{\Omega }_i}}\left[\left(\frac{{\left(\eta \beta \right)}_{xx}}{2}-{\eta }_x{\beta }_x\right)\left({\left|z\right|}^2+{\left|w\right|}^2\right)\right]\text{d}x\\ \mathrm{<mtext>\; \hspace{0.17em}\; </mtext>}\le \mathrm{<mtext>\; \hspace{0.17em}\; </mtext>}{C}_1{\Vert \phi \Vert }^2,\end{array}$

where $C_1=\underset{x}{\max }4\left\{\frac{{\left(\eta \beta \right)}_{xx}}{2}-{\eta }_x{\beta }_x\right\}$.

$I_2={\displaystyle {\int }_{{\Omega }_i}}\left[i\left({\bar{V}}_{PML}-V_{PML}\right)\bar{\phi }\phi \right]\text{d}x=2{\displaystyle {\int }_{{\Omega }_i}}\left[\mathrm{Im}\left\{V_{PML}\right\}\bar{\phi }\phi \right]\text{d}x\le C_2{\Vert \phi \Vert }^2,$

where $C_2={\eta }_x{\beta }_x-{\left(\eta \beta \right)}_{xx}$.

Applying the Lemmas 1 and 2, we have

${}_0{}^{C}D{}_t^{\alpha }{\Vert \phi \Vert }^2\le C_3{\Vert \phi \Vert }^2$

with $C_3=\underset{x}{\max }\left\{{\left(\eta \beta \right)}_{xx}-{\eta }_x{\beta }_x\right\}$, and

${\Vert \phi \left(\cdot \mathrm{,}t\right)\Vert }^2\le K\left(T\right){\Vert \phi \left(\cdot \mathrm{,0}\right)\Vert }^2$

in any interval $0\le t\le T$. This completes the proof. □

3. Finite Difference Scheme for Reduced IBVP

In this section, the following difference scheme is established for the equation and the following assumptions are made. The interval $\left[x_l-d\mathrm{,}{x}_r+d\right]$ is uniformly divided into M equal parts, the interval $\left[\mathrm{0,}T\right]$ is uniformly divided into N equal parts, the spatial step size is $h=\left(x_r-x_l+2d\right)/M$, the temporal step size is $\tau =T/N$, let $x_j=jh\left(0\le j\le M\right)$, $t_n=n\tau \left(0\le n\le N\right)$,

${\Omega }_h=\left\{x_j\mathrm{<mo>\; |\; </mo>0}\le j\le M\right\}$, ${\Omega }_{\tau }=\left\{t_n\mathrm{<mo>\; |\; </mo>0}\le n\le N\right\}$.

Assume $\Phi =\left\{{\varphi }_j^n\mathrm{<mo>\; |\; </mo>0}\le j\le M\mathrm{,0}\le n\le N\right\}$ is a grid function defined on ${\Omega }_h\times {\Omega }_{\tau }$, where ${\varphi }_j^n\approx \psi \left(x_j\mathrm{,}{t}_n\right)$. The following notations are reported

$\begin{array}{l}{\gamma }_k=\gamma \left(x_k\right),\mathrm{<mtext>\; \hspace{0.17em}\; </mtext>}{\delta }_x{\varphi }_j^n=\frac{{\varphi }_{j+\frac{1}{2}}^n-{\varphi }_{j-\frac{1}{2}}^n}{h},\mathrm{<mtext>\; \hspace{0.17em}\; </mtext>}{\varphi }_j^{n+\frac{1}{2}}=\frac{{\varphi }_j^{n+1}+{\varphi }_j^n}{2},\\ \mathrm{<mtext>\; \hspace{0.17em}\; </mtext>}{D}_{\tau }{\varphi }_j^n=\frac{{\tau }^{-\alpha }}{\Gamma \left(2-\alpha \right)}\left[a_0^{\left(\alpha \right)}{\varphi }_j^n-{\displaystyle \underset{l=1}{\overset{n-1}{\sum }}}\left(a_{n-l-1}^{\left(\alpha \right)}-a_{n-l}^{\left(\alpha \right)}\right){\varphi }_j^l-a_{n-1}^{\left(\alpha \right)}{\varphi }_j^0\right],\end{array}$

where $a_l^{\left(\alpha \right)}={\left(l+1\right)}^{1-\alpha }-l^{1-\alpha }$.

The PML function (9) can have the following approximation by applying the Crank-Nicolson finite difference scheme

$i{D}_{\tau }{\varphi }_j^n+{\gamma }_j{\delta }_x\left({\gamma }_j{\delta }_x{\varphi }_j^{n-\frac{1}{2}}\right)=\mathrm{0,}$ (19)

${\varphi }_j^0={\psi }_0\left(x_j\right),$ (20)

where $0\le j\le M$, and $0\le n\le N$.

4. Numerical Results

In this section, the absorbing function is given $\gamma \left(x\right)={\sigma }_0{\left(x-d\right)}^2$, where ${\sigma }_0>0$ is the absorption strength factor. Some numerical results are given to demonstrate the effectiveness of the PML approach in this section.

Example 1. Consider the time fractional Schrödinger equation with the initial condition

$\psi \left(x\mathrm{,0}\right)=\exp \left(-x^2+ikx\right)\mathrm{,}$

where k is the wave number. The computational interval is $\left[x_l-d\mathrm{,}{x}_r+d\right]=$

Figure 1. The influence of the thickness of the PML absorbing layers.

$\left[-5-d\mathrm{,5}+d\right]$, the wave number is chosen as $k=20$. Since the exact solution is undiscovered, we calculate the TFSE on a very fine grid and in a large area $\left[-\mathrm{15,15}\right]$ as reference solution.

Figure 1 plots the reference solution and the numerical solution with different PML width $d=0.2,0.4,0.7,1$. One can see that the performance of PML is greatly improved when the absorbing layer is enlarged. The width of PML is selected as $d=1$ in the following computations.

Figure 2 shows the numerical and exact solutions for different fractional orders $\alpha =0.4,\mathrm{<mtext>\; \hspace{0.17em}\; </mtext>\; 0.8}$ at different times, and one can observed that the numerical solutions fit the exact solutions very well, which means that the designed absorbing function in PML is very efficient. Figure 3 studies the influence of the absorption strength factor, we can see the PML approach present approximate solution very well as ${\sigma }_0=1$ in this example, which demonstrates the validity of the proposed method.

Figure 2. Numerical solutions compared with exact solutions for different orders $\alpha =0.4$ and $\alpha =0.8$.

(a)(b)

Figure 3. The numerical solutions with different ${\sigma }_0$ at different times.

5. Conclusion

The perfectly matched layer function of the time fractional Schrödinger equation on an unbounded domain is developed by applying the perfectly matched layer approach to obtain an initial boundary value problem on a bounded computational domain, which can be solved efficiently by adopting finite difference method. Based on the idea of smooth exterior scaling, the stability of the reduced initial boundary value problem is presented. The numerical results are reported to demonstrate the validity of the proposed method.

Acknowledgements

This research is supported by Shandong provincial natural science foundation (Grant No. ZR2019BA002).

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

References

[1]	Feynman, R. and Hibbs, A. (1965) Quantum Mechanics and Path Integrals. McGraw Hill, New York.
[2]	Laskin, N. (2000) Fractional Quantum Mechanics. Physical Review E, 62, 3135-3145. https://doi.org/10.1103/PhysRevE.62.3135
[3]	Laskin, N. (2000) Fractional Quantum Mechanics and Lévy Path Integrals. Physics Letters A, 268, 298-305. https://doi.org/10.1016/S0375-9601(00)00201-2
[4]	Naber, M. (2004) Time Fractional Schrödinger Equation. Journal of Mathematical Physics, 45, 3339-3352. https://doi.org/10.1063/1.1769611
[5]	Sjögreen, B. and Petersson, N.A. (2005) Perfectly Matched Layers for Maxwell’s Equations in Second Order Wave Equations. Journal of Computational Physics, 209, 19-46. https://doi.org/10.1016/j.jcp.2005.03.011
[6]	Singer, I. and Turkel, E. (2004) A Perfectly Matched Layer for the Helmholtz Equation in a Semi-Infinite Strip. Journal of Computational Physics, 201, 439-465. https://doi.org/10.1016/j.jcp.2004.06.010
[7]	Karim, D., Ballandras, S., Laroche, T., et al. (2013) Finite Element Analysis in Combination with Perfectly Matched Layer to the Numerical Modeling of Acoustic Devices in Piezoelectric. Applied Mathematics, 4, 64-71. https://doi.org/10.4236/am.2013.45A008
[8]	Kusnezov, D., Bulgac, A. and Dang, G. (1999) Quantum Lévy Processes and Fractional Kinetics. Physical Review Letters, 82, 1136. https://doi.org/10.1103/PhysRevLett.82.1136
[9]	Lomin, A. (2000) Fractional-Time Quantum Dynamics. Physical Review E, 62, 3135-3145. https://doi.org/10.1103/PhysRevE.62.3135
[10]	Wang, S. and Xu, M. (2007) Generalized Fractional Schrödinger Equation with Space-Time Fractional Derivatives. Journal of Mathematical Physics, 48, Article ID: 043502. https://doi.org/10.1063/1.2716203
[11]	Tarasov, V. (2006) Fractional Heisenberg Equation. Physics Letters A, 372, 2984-2988. https://doi.org/10.1016/j.physleta.2008.01.037
[12]	Duo, S. and Zhang, Y. (2015) Computing the Ground and First Excited States of Fractional Schrödinger Equations in an Infinite Potential Well. Communications in Computational Physics, 18, 321-350. https://doi.org/10.4208/cicp.300414.120215a
[13]	Zheng, X. and Wang, H. (2022) Analysis and Discretization of a Variable-Order Fractional Wave Equation. Communications in Nonlinear Science and Numerical Simulation, 104, Article ID: 106047. https://doi.org/10.1016/j.cnsns.2021.106047
[14]	Khan, N., Jamil, M. and Ara, A. (2012) Approximate Solutions to Time-Fractional Schrödinger Equation via Homotopy Analysis Method. ISRN Mathematical Physics, 2012, Article ID: 197068. https://doi.org/10.5402/2012/197068
[15]	Odibat, Z., Momani, S. and Alawneh, A. (2008) Analytic Study on Time-Fractional Schrödinger Equations: Exact Solutions by GDTM. Journal of Physics: Conference Series, 96, Article ID: 012066. https://doi.org/10.1088/1742-6596/96/1/012066
[16]	Berenger, J.P. (1994) A Perfectly Matched Layer for the Absorption of Electromagnetic Waves. Journal of Computational Physics, 114, 185-200. https://doi.org/10.1006/jcph.1994.1159
[17]	Narayanan, T. and Yoder, P. (2007) Perfectly Matched Layer Boundary Conditions for Quantum Phase Space Transport. Physics Letters A, 367, 288-290. https://doi.org/10.1016/j.physleta.2007.03.052
[18]	Cohen, G. and Imperiale, S. (2012) Perfectly Matched Layer with Mixed Spectral Elements for the Propagation of Linearized Water Waves. Communications in Computational Physics, 11, 285-302. https://doi.org/10.4208/cicp.201109.261110s
[19]	Araujo, E. and Pestana, R. (2020) Perfectly Matched Layer Boundary Conditions for the Second-Order Acoustic Wave Equation Solved by the Rapid Expansion Method. Geophysical Prospecting, 68, 572-590. https://doi.org/10.1111/1365-2478.12868
[20]	Zheng, C. (2007) A Perfectly Matched Layer Approach to the Nonlinear Schrödinger Wave Equations. Journal of Computational Physics, 227, 537-556. https://doi.org/10.1016/j.jcp.2007.08.004
[21]	Nissen, A. and Kreiss, G. (2011) An Optimized Perfectly Matched Layer for the Schrödinger Equation. Communications in Computational Physics, 9, 147-179. https://doi.org/10.4208/cicp.010909.010410a
[22]	Dohnal, T. (2009) Perfectly Matched Layers for Coupled Nonlinear Schröinger Equations with Mixed Derivatives. Journal of Computational Physics, 228, 8752-8765. https://doi.org/10.1016/j.jcp.2009.08.023
[23]	Antoine, X., Lorin, E. and Tang, Q. (2017) A Friendly Review of Absorbing Boundary Conditions and Perfectly Matched Layers for Classical and Relativistic Quantum Waves Equations. Molecular Physics, 115, 1861-1879. https://doi.org/10.1080/00268976.2017.1290834
[24]	Li, D. and Zhang, J. (2016) Efficient Implementation to Numerically Solve the Nonlinear Time Fractional Parabolic Problems on Unbounded Spatial Domain. Journal of Computational Physics, 322, 415-428. https://doi.org/10.1016/j.jcp.2016.06.046
[25]	Zhang, J., Li, D. and Antoine, X. (2019) Efficient Numerical Computation of Time-Fractional Nonlinear Schrödinger Equations in Unbounded Domain. Communications in Computational Physics, 25, 218-243. https://doi.org/10.4208/cicp.OA-2017-0195
[26]	Antoine, X. and Lorin, E. (2019) Towards Perfectly Matched Layers for Time-Dependent Space Fractional PDEs. Journal of Computational Physics, 391, 59-90. https://doi.org/10.1016/j.jcp.2019.04.025
[27]	Antoine, X., Lorin, E. and Zhang, Y. (2021) Derivation and Analysis of Computational Methods for Fractional Laplacian Equations with Absorbing Layers. Numerical Algorithms, 87, 409-444. https://doi.org/10.1016/j.jcp.2019.04.025
[28]	Li, D., Zhang, J. and Zhang, Z. (2018) The Numerical Computation of the Time Fractional Schrödinger Equation on Unbounded Domain. Computational Methods in Applied Mathematics, 18, 77-94. https://doi.org/10.1515/cmam-2017-0038
[29]	Podlubny, I. (1999) Fractional Differential Equations. Academic Press, San Diego.