2017-09-14 59 views
3

我想在python中使用scipy.odeint函數解決以下方程。在Python中使用odeint實現一個集成數學方程

equation 1

目前我可以使用下面的腳本來實現這種形式的方程

equation 2

在蟒:

def dY(y1, x): 
    a = 0.001 
    yin = 1 
    C = 0.01 
    N = 1 
    dC = C/N 
    b1 = 0 
    return (a/dC)*(yin-y1)+b1*dC 

x = np.linspace(0,20,1000) 
y0 = 0 
res = odeint(dY, y0, x) 
plt.plot(t,res, '-') 
plt.show() 

我與第一個方程問題'一世'。我不知道如何對方程進行積分,仍然能夠提供當前和以前的'y'(yi-1和yi)值。 '我'只是一個在0..100範圍內的序列號。

編輯1:

原始等式爲: Original equation

哪我改寫使用Y,X,A,B和C

EDIT2: 我編輯皮埃爾德Buyl」代碼和改變了N值。幸運的是我有一個驗證表來驗證結果。不幸的是,結果並不相同。

這裏是我的驗證表:

Validation image

這裏是numpy的輸出:

numpy output

使用的代碼:

def dY(y, x): 
    a = 0.001 
    yin = 1 
    C = 0.01 
    N = 3 
    dC = C/N 
    b1 = 0.01 
    y_diff = -np.copy(y) 
    y_diff[0] += yin 
    y_diff[1:] += y[:-1] 
    return (a/dC)*(y_diff)+b1*dC 

x = np.linspace(0,20,11) 
y0 = np.zeros(3) 
res = odeint(dY, y0, x) 
plt.plot(x,res, '-') 

正如你所看到的值與a不同n的偏移量爲0.02 ..

我是否錯過了導致此偏移的原因?

+0

這符號可能意味着你使用的歐拉前向整合爲一個單一的ODE爲Y或整合多個耦合常微分方程。哪一個? – duffymo

+0

@duffymo他們正在使用歐拉正向積分求解y –

+0

@MD'這個微分方程是否正確?你能否從你從哪裏得到一個參考? –

回答

1

的方程是「耦合的」常微分方程(見「常微分方程的系統」上Wikipedia

變量是含有y[0]y[1]等爲了解決ODE必須喂矢量作爲矢量。初始條件和功能dY必須返回一個向量以及

我已經修改您的代碼以實現這一結果:

def dY(y, x): 
    a = 0.001 
    yin = 1 
    C = 0.01 
    N = 1 
    dC = C/N 
    b1 = 0 
    y_diff = -np.copy(y) 
    y_diff[0] += yin 
    y_diff[1:] += y[:-1] 
    return (a/dC)*y_diff+b1*dC 

我有WRI將y[i-1] - y[i]部分作爲NumPy向量操作,並且特殊地使用座標y[0](即符號中的y1,但Python中的數組從0開始)。

的解決方案,使用的初始值爲0所有yi是

x = np.linspace(0,20,1000) 
y0 = np.zeros(4) 
res = odeint(dY, y0, x) 
plt.plot(x,res, '-') 
+0

謝謝你的實施。但是,我通過對y0和N進行了一些修改來執行代碼,以便輸出與我的驗證表的結果相匹配。請參閱Edit2,其中詳細說明了我的內容。 –