时间依赖的薛定谔方程奇偶解

问题描述 投票:1回答:1

我正在解决谐波势和初始高斯波函数下时间依赖的薛定谔方程的演化。处理hcut=12m=1,并在实部和虚部中分离波函数,根据实部和虚部分别得到两个耦合方程,分别称为yryc

xrange[xi,xf] trange[0,tf]

我使用的方法是:

首先在实部和虚部中分离波函数,即yr(x,t)和yc(x,t)。然后处理hcut = 2m = 1,并将波函数写为yr[x,t]+i*yc[x,t],我们从TDSE得到两个耦合方程。

1.D[yr[x,t],t]=-D[yc[x,t],x,x]+V[x]*yc[x,t]
2.D[yc[x,t],t]=D[yc[x,t],x,x]-V[x]*yc[x,t]

然后我将初始波函数指定为

yr[x,0]=exp[-x^2]
yc[x,0]=0

之后,使用有限差分格式,我试图从y [x,0]找到y [x,t],即,

yr[x,t+d]=yr[x,t]+d*D[yr[x,t],t]
        =yr[x,t]+d*(D[yc[x,t],x,x]+V[x]*yc[x,t])

在代码中索引为“a”,在复杂部分索引为“b”。

y[x[i],t[j]]的值存储为y[i,j],因为数组不能具有实数索引。

我使用的代码如下:

    function v(x) result(s)
    real::s,x
    s=x**2
    end function v

    real::t(10000),x(10000),yc(10000,10000),yr(10000,10000),tf,xi,xf,d
    integer::i,j,k,l,m,n,time
    write(*,*) "time of plot divided by step size"
    read(*,*) time

    tf=50
    xi=-10
    xf=10
    d=0.01
    x(1)=xi
    t(1)=0

    i=1
    do while(x(i).lt.xf)        !input all values of x in x(i) array
    x(i+1)=x(i)+d
    i=i+1
    end do

    n=1
    do while(t(n).lt.tf)        !input all values of t in t(i) array
    t(n+1)=t(n)+d
    n=n+1
    end do




    do j=1,i
    yr(j,1)=exp(-(x(j)-5)**2)       !input of initial wavefunction    
    yc(j,1)=0
    end do

    !calculation of wavefunction at higher time using finite element technique[y[x,t+d]=y[x]+d*D[y[x,t],t] and then replacing the partial derivative with time 
    !using equation 1 and 2 .

    l=1
    do while(l.le.i-2)
    k=1
         do while(t(k).lt.tf)
            yr(l,k+1)=yr(l,k)-(yc(l+2,k)-2*yc(l+1,k)+yc(l,k))/d&    
            +v(x(l))*yc(l,k)*d                                                                             

            yc(l,k+1)=yc(l,k)+(yr(l+2,k)-2*yr(l+1,k)+yr(l,k))/d&   
            -v(x(l))*yr(l,k)*d

            k=k+1

         end do
    l=l+1
    end do




    open(1,file="q.dat")
    do m=1,i-2
    write(1,*) x(m),(yr(m,time))**2+(yc(m,time))**2
    end do
    close(1)

    end

预期结果:y(x,t)^ 2 = yr(x,t)^ 2 + yc(x,t)^ 2

错误:波函数没有保持规律,只有t=0.050.06后,波函数变大,最大值变为e30的数量级,尽管高斯形状几乎保持不变,正如预期的那样,只有0.05秒已经过去了。

fortran physics gfortran
1个回答
-1
投票

这不是一个真正的编码问题,而是一个数学数学问题,我建议你先学习一些关于Schrödinger方程和类似偏微分方程的数值方法的教程,然后在专门讨论scientific computation的网站上提出更多要求。

首先,您选择的方法不能保留解决方案的规范。这可以通过每个时间步骤的标准化来保存,但它非常难看。

其次,您选择的方法是显式的(实际上,它可能是无条件不稳定的前向Euler方法)。这意味着可延迟的时间步骤将受到类似扩散的术语的严重限制(尽管这里很复杂)。一个简单的隐式方法也很好地保存了标准是Crank-Nicolson method。这是一种隐式方法,您必须在每个时间步长求解一个线性方程组。在1D中它非常简单,因为系统是三对角的。还有一些explicit schemes that may work,但他们更多参与,而不是你尝试过的天真计划。

第三,您不必将系统拆分为实部和虚部,可以使用复数来完成。

第四,您应该将空间网格步长dx和时间步长dt保持为具有已知值的单独变量。你可以计算最终系数d作为它们的比率,但你应该至少知道你的值od dx和dt。

第五,您不必在所有时间步骤中存储解决方案的所有值。在较大问题所需的记忆时间内,很快就会成为禁令。存储最后一个时间步和当前时间步就足够了。更多用于多步骤方法和Runge-Kutta方法的一些辅助步骤。

© www.soinside.com 2019 - 2024. All rights reserved.