0
我想重現下面顯示的兩個圖,它們是作爲時間函數的溫度和濃度分佈圖。我檢查了我的方法和代碼一百萬次,但似乎無法找到它的錯誤,但我無法再現這些圖。除了CA和T之外,所有的值都是常量。它可能是來自scipy的odeint的準確性問題嗎?任何幫助將非常感激!用於ODE系統的scipy odeint的錯誤
這兩個方程如下:
DCA/DT = Q *(CAI - CA)/ V - K * CA
的dT/dt的= W *(鈦 - T)/(V p)+ d_HRķCA /(p C)+ UA *(TC - T)/(V p C)
的代碼:
import numpy as np
from scipy.integrate import odeint
import matplotlib.pyplot as plt
def ODESolve(y, t, q, CAi, V, k0, w, Ti, p, dH_R, C, UA, Tc):
CA, T = y
k = k0*np.exp(8750*1/T)
dydt = [q*(CAi - CA)/V - k*CA, w*(Ti - T)/(V*p) + \
dH_R*k*CA/(p*C) + UA*(Tc - T)/(V*p*C)]
return dydt
q = 100
CAi = 1.0
V = 100
p = 1000
C = .239
dH_R = 5*(10**4)
k0 = 7.2*(10**10)
UA = 5*10**4
CA0 = .5
T0 = 350
Ti = T0
w = p*q
y0 = [CA0, T0]
t = np.linspace(0, 20, 100)
Tc = 305
sol1 = odeint(ODESolve, y0, t, args = (q, CAi, V, k0, w, Ti, p, dH_R, C, UA, Tc))
Tc = 300
sol2 = odeint(ODESolve, y0, t, args = (q, CAi, V, k0, w, Ti, p, dH_R, C, UA, Tc))
Tc = 290
sol3 = odeint(ODESolve, y0, t, args = (q, CAi, V, k0, w, Ti, p, dH_R, C, UA, Tc))
plt.figure(1)
plt.plot(t, sol1[:,0], label = 'Tc = 305')
plt.plot(t, sol2[:,0], label = 'Tc = 300')
plt.plot(t, sol3[:,0], label = 'Tc = 290')
plt.ylim(ymax = 1, ymin = 0)
plt.title ('CA(t)')
plt.legend(loc = 'best')
plt.figure(2)
plt.plot(t, sol1[:,1], label = 'Tc = 305')
plt.plot(t, sol2[:,1], label = 'Tc = 300')
plt.plot(t, sol3[:,1], label = 'Tc = 290')
plt.ylim(ymax = 450, ymin = 300)
plt.legend(loc = 'best')
plt.title ('T(t)')
plt.show()
這裏是圖應該生產什麼:
這裏是我的代碼的輸出上面:
仔細檢查此公式:'k = k0 * np.exp(8750 * 1/T)'。在一個完整的黑暗中,我改變了指數內的符號,所以'k = k0 * np.exp(-8750 * 1/T)',得到的結果看起來類似於你預期的情節顯示。 –
你是上帝......這是完全正確的!哇,真是一種解脫。如果您將此作爲答案提交,我很樂意提交您的答案作爲答案。謝謝! –
好的,提交答案。 –