from pylab import *
figure(3)
n = 39
psi0= zeros(n) #
psi = psi0 #
x = linspace(0,1,n)
V = zeros(n)
#V = 39*pow(x,2)
for i in arange(n):
if i<n/2.:
V[i] = 73.
else:
V[i] = 0.
psi0[0] = 1.
psi0[1] = psi0[0] #for odd function, use another value,
#but psi0[0] must be zero
t = 0
dx = 1./n
E = 1.
dE = .1
err = .05
while t< 771:
#k = 2*dx*dx*(E-V)
for i in range (1,n-1):
k = 2*dx*dx*(E-V[i])
psi[i+1] = 2*psi0[i]-psi0[i-1]-k*psi0[i]
psi0 = psi
if abs(psi[n-1])<err:
print E
plot(x,psi)
t += 1
E += dE
xlabel('x')
ylabel('psi')
title('Fungsi Gelombang')
grid(True)
savefig("shooting.png")
figure(2)
plot(x,V)
xlabel('x')
ylabel('V')
title('Potensial')
grid(True)
show()
Here's some result with different potential V
No comments:
Post a Comment