并不是太理解您的意思，但我还是让电脑算了一下，结果和以前一样（昨天可能运行其他内容，所以速度较慢吧，今天挺快的）。库朗条件貌似没有应用到复数范围内的吧？fdtd倒是有篇论文，或许可以借鉴一下，正在研究

你是用迭代法求解，你只迭代了一次，对于AX=B这样的方程至少要迭代A的阶数次才能得到近似解，你的那两层循环只是把边值在整个区域“匀和”了一次，你需要“匀和”N次结果才能比较接近解
或说回来，如果你先形成细数矩阵A，然后用直接法求解，X=A^{-1}B就不会有问题

kerb 发表于 2017-1-5 22:39
do j=1,n-1
        do i=2,m-1
            x(i)=i*h

加个k的目的是什么？感觉除了一直算下去、然后输出同样的结果，并没有什么区别呀？另外，电脑已经算10分钟了
明天起床看看库朗条件吧，谢谢啦

[Fortran] 纯文本查看 复制代码

program time_dependent_schrodinger_equation
    implicit none
    integer,parameter::M=150,N=500
    integer::i,j
    real*16::hba=1.,mu=1./2
!hba即（普朗克常数/2/pi），mu为质量
    real*16::l=150,T1=5
    real,external::V ！势能函数
    real*16,dimension(M)::x,x1
    real*16,dimension(N)::t,normalize,psy4
    real*16,dimension(M,N)::psy2,psy3  ！psy2为对复数psy取模并求其平方的大小，psy3，即为归一化后的波恩诠释下的概率。
    complex*16,dimension(M,N)::psy
    complex::cj
    real*16::h,tao,pi
    h=l/M
    tao=T1/N
    pi=asin(1.)*2
    cj=(0.,1.)  !虚数单位
!给出初始边界及初始值
       do j=1,N
        psy(1,j)=0
        psy(M,j)=0 ！边界设定为无限深势阱
    end do
    do i=2,M-1
        x1(i)=i*h
        psy(i,1)=exp(((-1)*(x1(i)*0.2-5)**2./2.)/(2.*pi))   ！阱内给出初始值
    end do                          !x(i)比着正常多乘上个0。2，平移了波的位置
!中心差法解薛定谔方程
    do j=1,N-1
        do i=2,m-1
            psy(i,j+1)=psy(i,j)+cj*tao*((psy(i+1,j)-2*psy(i,j)+psy(i-1,j))/(h**2)-V(i)*psy(i,j))    ！核心内容，将中心差法代入薛定谔方程并整理
            psy2(i,j+1)=(abs(psy(i,j+1)))**2
       !     write(13,*) psy2(i,j+1)
        end do
    end do
!归一化
    do j=2,n
        normalize(j)=(sum(psy2(1:m,j)))!**0.5   
    end do
    do j=2,N
        do i=2,m-1
            x(i)=i*h
            t(j)=j*tao
            psy3(i,j)=psy2(i,j)/(normalize(j))!**2
            write(13,*) x(i),t(j),psy3(i,j)
        end do
        psy4(j)=sum(psy3(1:M,j))
    !    write(14,*) psy4(j)   !检验归一化的正确性
    end do
    stop
    end
!给出势函数
real    function V(i)
        implicit none
        integer::i
        v=0
        if (abs(i-75)-5<1.e-6) then
          V=10.
        end if
!   if abs(x-40)<3&abs(x-60)<3
!   V=10
        return
end

这样求解不行
试一下这样

[Fortran] 纯文本查看 复制代码

do j=1,n-1
        do i=2,m-1
            x(i)=i*h
            t(j+1)=(j+1)*tao
            y(i,j+1)=y(i,j)+tao*((y(i+1,j)-2*y(i,j)+y(i-1,j))/(h**2))
            write(22,*) x(i),t(j+1),y(i,j+1)
        end do
    end do

你搜一下courant condition
双曲类偏微分方程差分法求解，对dx,dy,dt等有一定要求，你自己搜一下

会不会边界条件已经暗示了我的波函数只能是驻波？
不知道为何，就连我给的热传导程序，也都是个静止的图像。

[Fortran] 纯文本查看 复制代码

do k=1,m*n
do j=1,n-1
        do i=2,m-1
            x(i)=i*h
            t(j+1)=(j+1)*tao
            y(i,j+1)=y(i,j)+tao*((y(i+1,j)-2*y(i,j)+y(i-1,j))/(h**2))
            if(k==m*n)write(22,*) x(i),t(j+1),y(i,j+1)
        end do
    end do
enddo

kerb 发表于 2017-1-3 20:11
时间步长选着不合适
你可以找本时域有限差分法（FDTD）看看，
附件是一个解释

我根据c（光速）=dx/dt作了图。M=15000,N=50，l=1500,T1=0.05，可是图像也并未移动。另外，在我的程序里，我应该是给光速暗含了重定义，并不等与传统的10**8量级，所以应该并不作影响吧？

时间步长选着不合适
你可以找本时域有限差分法（FDTD）看看，
附件是一个解释