2016-11-10 47 views
2

我需要改变我的绘图在薛定谔方程,y轴上的比例,以显示理论计算与我们的差异,大约为0.01%的差异。所以在剧情上我得到的规模还不够小,以显示差异。这是我的项目代码。薛定谔方程的变比缩放

# -*- coding: utf-8 -*- 
""" 
Created on Sat Nov 05 12:25:14 2016 

@author: produce 
""" 
from __future__ import print_function 
import numpy as np 
import matplotlib.pyplot as plt 

# 
c = .5/500 # c = delta x 
x = np.arange(0, .5, c) # creates array of argument values from 0 to 1/2 in increments 
#      of delta x = c 
psi = np.zeros(len(x)) # creates array of zeros which will be replaced by y values 

k = 20 # starting energy for calculator of E 
ans = 0 # The value of k, when we have y as between 0.004 and 0 
ansPsi = 0 
diff = 0.001 
increment = 0.0001 
done = False 
while 1: 
    # print k 
    psi[0] = 1 
    psi[1] = 1 
    for i in range(0, len(x) - 2): 
     psi[i + 2] = psi[i + 1] + (psi[i + 1] - psi[i]) - 2 * k * c * c * psi[i] 
    # plt.plot(x,psi) 
    # print(x,psi) 
    # print (psi[i+2]--->) 
    if (float(psi[i + 2]) < 0.004 and float(psi[i + 2]) > 0): 
     ans = k 
     ansPsi = psi[i + 2] 
     # print ("NOW ENTERING INNER LOOP") 
     while 1: # would be an infinite loop, but have a break statement 
      # k = k - 0.00001 
      k = k + increment 
      for i in range(0, len(x) - 2): 
       psi[i + 2] = psi[i + 1] + (psi[i + 1] - psi[i]) - 2 * k * c * c * psi[i] 
      plt.plot(x, psi, 'r') #red solid line 
      if (psi[i + 2] > ansPsi or psi[i + 2] < 0): 
       done = True 
       break 
      else: 
       ansPsi = psi[i + 2] 
       ans = k 
       # print (k, psi[i+2]) 

    if done: 
     break 
    k = k - diff 

print("Value of k:", ans, "Value of Y:", ansPsi) # prints our answer for energy and psi[1/2] 
k1 = 10 # 1st Higher Energy Value 
k2 = 7 # 2nd Higher Energy Value 
k3 = 3 # 1st Lower Energy Value 
k4 = 1 # 2nd Lower Energy Value 
kt = np.pi * np.pi * .5 # theoretical value 

psi1 = np.zeros(len(x)) 
psi1[0] = 1 
psi1[1] = 1 
for i in range(0, len(x) - 2): 
    psi1[i + 2] = psi1[i + 1] + (psi1[i + 1] - psi1[i]) - 2 * k1 * c * c * psi1[i] 


# psi2 = np.zeros(len(x)) 
# psi2[0] = 1 
# psi2[1] = 1 
# for i in range (0,len(x)-2): 
# psi2[i+2] = psi2[i+1] + (psi2[i+1] - psi2[i]) - 2*k2*c*c*psi2[i] 
# plt.plot(x,psi2,'k') 

# psi3 = np.zeros(len(x)) 
# psi3[0] = 1 
# psi3[1] = 1 
# for i in range (0,len(x)-2): 
# psi3[i+2] = psi3[i+1] + (psi3[i+1] - psi3[i]) - 2*k3*c*c*psi3[i] 
# plt.plot(x,psi3,'p') 

psi4 = np.zeros(len(x)) 
psi4[0] = 1 
psi4[1] = 1 
for i in range(0, len(x) - 2): 
    psi4[i + 2] = psi4[i + 1] + (psi4[i + 1] - psi4[i]) - 2 * k4 * c * c * psi4[i] 

plt.plot(x, psi, 'r-', label='Corrected Energy') 

psiT = np.zeros(len(x)) 
psiT[0] = 1 
psiT[1] = 1 
for i in range(0, len(x) - 2): 
    psiT[i + 2] = psiT[i + 1] + (psiT[i + 1] - psiT[i]) - 2 * kt * c * c * psiT[i] 
plt.plot(x, psiT, 'b-', label='Theoretical Energy') 
plt.ylabel("Value of Psi") 
plt.xlabel("X value from 0 to 0.5") 
plt.title("Schrodingers equation for varying inital energy") 
plt.legend(loc=3) 
plt.yscale() 
plt.show() 

回答

0

您共享的代码由于plt.yscale()需要参数而失败。我只是评论说,线。

由于您的理论能量曲线和修正后的能量曲线相差很小,所以不可能缩放y轴仍然可以看到x的整个范围内的两条曲线(即 - 从0到0.5) 。相反,也许你应该绘制两条曲线的区别?

plt.plot(x, psiT-psi) 
plt.title("Size of Correction for Varying Initial Energy") 
plt.ylabel(r"$\Delta$E") 
plt.xlabel("X value from 0 to 0.5") 
plt.show() 

另外,在x和y标签上贴一些单位可能会很好。 :)