二階常微分方程RK4

我需要繪製顯示2個變量的圖形，具有二階常微分方程RK4，到目前爲止，我已經做到了這一點二階常微分方程RK4

from numpy import arange 
from pylab import plot,xlabel,ylabel,show 
Qger = 400 
K = 20 
T1 = 150 
T2 = 60 
N = 1000 
h = (T2-T1)/N 
rpoints = arange(6.0,8.0,h) 
xpoints = [] 
x = 423 
def df(s,t): 
    dTdt = -Qger*t/(2*K) + 172.8/t 
    return dTdt 

for r in rpoints: 
    xpoints.append(x) 
    k1 = h*df(x,r) 
    k2 = h*df(x+0.5*k1,r+0.5*h) 
    k3 = h*df(x+0.5*k2,r+0.5*h) 
    k4 = h*df(x+k3,r+h) 
    x += (k1+2*k2+2*k3+k4)/6 
pylab.plot(rpoints,xpoints) 
pylab.xlabel("Raio") 
pylab.ylabel("Temperatura") 
pylab.show

但是，這是一個RK4一階微分方程，因爲我不知道和集成手，但我不能這樣做，也不使用scipy，所以任何人都可以向我解釋如何集成此功能或使用RK4與二階ODE。該功能如下。

這是函數，只有T和R是變量，其餘爲0

來源

2017-06-21 Pedro Gabriel Carvalho

[odeint'的scipy文檔]（https://docs.scipy.org/doc/scipy-0.19.0/reference/generated/scipy.integrate.odeint.html）有一個解決方案的詳細示例的二階ODE。研究這個例子並回來（也許在另一個問題中......）以獲得關於特定問題的幫助。因爲你的問題是一個壞的... – gboffi

你應該能夠把上面的「半獨立」的形式，那是說的dT/dt僅就偏導數r而言。如果可以找到相當於dT/dt的數值或近似值，即dT/dt = df（r，...）的RHS，則顯式RK4可以適用。

在這種方法中，您的時間步進方法（RK4）僅適用於溫度相對於時間的一階導數。

來源

2017-07-11 23:14:01

二階常微分方程RK4

回答

相關問題