Python：odeint解決可變係數（QHO）ODE

Question

我試圖通過使用odeint數值求解方程y'' + (epsilon-x^2)y = 0。我知道解決方案（QHO的wavefunctions），但odeint的輸出與它沒有明顯的關係。我可以用常係數來解決ODE問題，但只要我選擇了可​​變參數，我就無法解決任何我嘗試過的問題。這裏是我的代碼：

#!/usr/bin/python2 
from __future__ import division 
import numpy as np 
import matplotlib.pyplot as plt 
import scipy.integrate as spi 

x = np.linspace(-5,5,1e4) 

n = 0 
epsilon = 2*n+1 

def D(Y,x): 
    return np.array([Y[1], (epsilon-x**2)*Y[0]]) 

Y0 = [0,1] 

Y = spi.odeint(D,Y0,x) 
# Y is an array with the first column being y(x) and the second y'(x) for all x 
plt.plot(x,Y[:,0],label='num') 
#plt.plot(x,Y[:,1],label='numderiv') 

plt.legend() 
plt.show() 

和劇情： [沒有足夠的代表：] https://drive.google.com/file/d/0B6840LH2NhNpdUVucUxzUGFpZUk/edit?usp=sharing

在這裏尋找的解決方案圖：http://hyperphysics.phy-astr.gsu.edu/hbase/quantum/hosc5.html

Answer 1

它看起來像你的公式不正確的解釋。你有一個微分方程y'' + (epsilon-x^2)y = 0，但你忘記了矢量形式的減號。特別是它應該是

y[0]' = y[1] 
y[1]' = -(epsilon-x^2)y[0] 

所以（添加減號在小量長期的正面

def D(Y,x): 
    return np.array([Y[1], -(epsilon-x**2)*Y[0]]) 

事實上，你有情節與DE y'' + (epsilon-x^2)y = 0一致看吧：Wolphram Alpha