四阶龙格库塔法解含时薛定谔方程,一维线性谐振子附加含时驱动力。

external f
parameter(n=100001,m=11)
Dimension Y(m),D(m),z(m,n),b(m),c(m,n)
complex(8) y,d,z,b,h
real(8) t,a
t=0.0d0

do j=1,m
Y(j)=(0.0,0.0)
enddo
Y(1)=(1.0d0,0.0d0)

h=2500.d0/(n-1)
call grkt1(t,y,m,h,n,z,f,d,b)
write(*,*)
call Psi(z)
end

subroutine Psi(z)
parameter(n=100001,m=11,pi=3.1415926d0)
complex(8) z(m,n),p(20,1000)
real(8) a(m,n)
integer t

open(200,File='test.txt')
open(201,File='test1.txt')
open(202,File='test2.txt')
open(203,File='test3.txt')

do k=1,20,1
do j=1,1000
x=-5.0d0+j*0.01d0
P(k,j)=(0.d0,0.d0)
do i=1,m
p(k,j)=p(k,j)+z(i,k*5000)/sqrt((sqrt(pi)*(2**i)*JC(i)))*exp(-x*x/2.0d0)*H(i,x*1.0)*exp((0.d0,-1.d0)*(2.d0*i+1.d0)*k*5000*(2500.d0/(n-1)))
!p(k,j)=p(k,j)+z(i,k)/sqrt((sqrt(pi)*(2**i)*JC(i)))*exp(-x*x/2.0d0)*H(i,x*1.0)*exp((0.d0,-1.d0)*(2.d0*i+1.d0)*k*(2500.d0/(n-1)))
enddo
enddo
enddo

do j=1,1000
write(200,*)(-5.d0+j*0.01d0),(abs(p(i,j))**2,i=1,2)
write(201,*)(abs(p(i,j))**2,i=3,5)
write(202,*)(abs(p(i,j))**2,i=6,8)
write(203,*)(abs(p(i,j))**2,i=9,11)
enddo

end subroutine

subroutine F(t,y,m,d)

Dimension Y(m),D(m)
complex(8) y,d
real(8) t
D(1)=-(0.d0,1.d0)*g(t)*(Y(2)*sqrt(1.d0/2.d0)*exp(-2.d0*(0,1)*t))
do j=2,m-1
D(j)=-(0,1)*g(t)*(Y(j-1)*sqrt((j-1)/2.d0)*exp(2.d0*(0,1)*t)+Y(j+1)*sqrt(j/2.)*exp(-2.d0*(0,1)*t))
enddo
D(m)=-(0,1)*g(t)*(Y(m-1)*sqrt((m-1)/2.d0)*exp(-2.d0*(0,1)*t))
return
end

function g(t)
real(8) g,t
pi=3.1415926535d0
if(t<=2000.d0)then
g=0.5*(1+sin(t/2000.*pi-pi/2.))
!    g=t/2000.d0
else
g=1.00d0
endif
end

real recursive function H(i,x) result(Hi)
select case(i)
case(1)
Hi=1.0
case(2)
Hi=2.0*x
case(3:)
Hi=2*x*h(i-1,x)-2*(i-1)*H(i-2,x)
end select
end

integer Recursive function JC(i) result(J)
if(i==0)then
J=1
else
J=i*JC(i-1)
endif
end

SUBROUTINE GRKT1(T,Y,M,H,N,Z,F,D,B)
DIMENSION Y(M),D(M),Z(M,N),A(4),B(M)
complex(8) Y,D,Z,A,B,H,X,TT
real(8) t
A(1)=H/2.0d0
A(2)=A(1)
A(3)=H
A(4)=H
DO I=1,M
Z(I,1)=Y(I)
enddo
X=T
DO J=2,N
CALL F(T,Y,M,D)
DO I=1,M
B(I)=Y(I)
enddo
DO K=1,3
DO I=1,M
Y(I)=Z(I,J-1)+A(K)*D(I)
B(I)=B(I)+A(K+1)*D(I)/3.0
enddo
TT=T+A(K)
CALL F(TT,Y,M,D)
enddo
DO I=1,M
Y(I)=B(I)+H*D(I)/6.0
enddo
DO I=1,M
Z(I,J)=Y(I)
enddo
T=T+H
enddo
T=X
RETURN
END

作者 hsyyf

《四阶龙格库塔法解含时薛定谔方程》有4条评论

回复 tusooa 取消回复

您的电子邮箱地址不会被公开。 必填项已用 * 标注